From 2e2e97cf30ff60d5ae941cc9229247c8683744ee Mon Sep 17 00:00:00 2001 From: David A Minton Date: Thu, 20 Apr 2023 13:24:55 -0400 Subject: [PATCH] Major restructing to move a bunch of basic utillity functions that operate on standard data types (DP, I4B, etc) to the base_module. This will make it easier to use them inside the coarray module. Also built coarray submodules to put implementations in. Restructuring still needs more thorough testing --- src/CMakeLists.txt | 2 + src/base/base_module.f90 | 1729 +++++++++++++++++++- src/coarray/coarray_clone.f90 | 494 ++++++ src/coarray/coarray_collect.f90 | 296 ++++ src/coarray/coarray_module.f90 | 992 ++---------- src/collision/collision_util.f90 | 2 +- src/encounter/encounter_check.f90 | 20 +- src/encounter/encounter_module.f90 | 6 +- src/encounter/encounter_util.f90 | 65 +- src/fraggle/fraggle_util.f90 | 4 +- src/rmvs/rmvs_coarray.f90 | 48 +- src/rmvs/rmvs_module.f90 | 17 + src/rmvs/rmvs_util.f90 | 104 +- src/swiftest/swiftest_coarray.f90 | 124 +- src/swiftest/swiftest_driver.f90 | 46 +- src/swiftest/swiftest_io.f90 | 6 +- src/swiftest/swiftest_module.f90 | 297 +--- src/swiftest/swiftest_util.f90 | 2370 +++++----------------------- src/symba/symba_util.f90 | 64 +- src/whm/whm_util.f90 | 62 +- 20 files changed, 3389 insertions(+), 3359 deletions(-) create mode 100644 src/coarray/coarray_clone.f90 create mode 100644 src/coarray/coarray_collect.f90 diff --git a/src/CMakeLists.txt b/src/CMakeLists.txt index 4e1caebba..d2a04991c 100644 --- a/src/CMakeLists.txt +++ b/src/CMakeLists.txt @@ -85,6 +85,8 @@ SET(FAST_MATH_FILES SET(COARRAY_FILES ${SRC}/coarray/coarray_module.f90 + ${SRC}/coarray/coarray_clone.f90 + ${SRC}/coarray/coarray_collect.f90 ${SRC}/swiftest/swiftest_coarray.f90 ${SRC}/whm/whm_coarray.f90 ${SRC}/rmvs/rmvs_coarray.f90 diff --git a/src/base/base_module.f90 b/src/base/base_module.f90 index 3c33a47c1..84e5a1f2e 100644 --- a/src/base/base_module.f90 +++ b/src/base/base_module.f90 @@ -220,8 +220,285 @@ end subroutine abstract_util_dealloc_object type, abstract :: base_nbody_system end type base_nbody_system + + interface util_append + module procedure base_util_append_arr_char_string + module procedure base_util_append_arr_DP + module procedure base_util_append_arr_DPvec + module procedure base_util_append_arr_I4B + module procedure base_util_append_arr_logical + end interface + + interface util_fill + module procedure base_util_fill_arr_char_string + module procedure base_util_fill_arr_DP + module procedure base_util_fill_arr_DPvec + module procedure base_util_fill_arr_I4B + module procedure base_util_fill_arr_logical + end interface + + interface util_resize + module procedure base_util_resize_arr_char_string + module procedure base_util_resize_arr_DP + module procedure base_util_resize_arr_DPvec + module procedure base_util_resize_arr_I4B + module procedure base_util_resize_arr_logical + end interface + + interface util_sort + module procedure base_util_sort_i4b + module procedure base_util_sort_index_i4b + module procedure base_util_sort_index_I4B_I8Bind + module procedure base_util_sort_index_I8B_I8Bind + module procedure base_util_sort_sp + module procedure base_util_sort_index_sp + module procedure base_util_sort_dp + module procedure base_util_sort_index_dp + end interface + + interface util_sort_rearrange + module procedure base_util_sort_rearrange_arr_char_string + module procedure base_util_sort_rearrange_arr_DP + module procedure base_util_sort_rearrange_arr_DPvec + module procedure base_util_sort_rearrange_arr_I4B + module procedure base_util_sort_rearrange_arr_I4B_I8Bind + module procedure base_util_sort_rearrange_arr_logical + module procedure base_util_sort_rearrange_arr_logical_I8Bind + end interface + + interface util_spill + module procedure base_util_spill_arr_char_string + module procedure base_util_spill_arr_DP + module procedure base_util_spill_arr_DPvec + module procedure base_util_spill_arr_I4B + module procedure base_util_spill_arr_I8B + module procedure base_util_spill_arr_logical + end interface + + interface util_unique + module procedure base_util_unique_DP + module procedure base_util_unique_I4B + end interface + contains + subroutine base_util_append_arr_char_string(arr, source, nold, lsource_mask) + !! author: David A. Minton + !! + !! Append a single array of character string type onto another. If the destination array is not allocated, or is not big enough, this will allocate space for it. + implicit none + ! Arguments + character(len=STRMAX), dimension(:), allocatable, intent(inout) :: arr !! Destination array + character(len=STRMAX), dimension(:), intent(in) :: source !! Array to append + integer(I4B), intent(in), optional :: nold !! Extent of original array. If passed, the source array will begin at arr(nold+1). Otherwise, the size of arr will be used. + logical, dimension(:), intent(in), optional :: lsource_mask !! Logical mask indicating which elements to append to + ! Internals + integer(I4B) :: nnew, nsrc, nend_orig + + if (present(lsource_mask)) then + nsrc = count(lsource_mask(:)) + else + nsrc = size(source) + end if + if (nsrc == 0) return + + if (.not.allocated(arr)) then + nend_orig = 0 + allocate(arr(nsrc)) + else + if (present(nold)) then + nend_orig = nold + else + nend_orig = size(arr) + end if + call util_resize(arr, nend_orig + nsrc) + end if + nnew = nend_orig + nsrc + + if (present(lsource_mask)) then + arr(nold + 1:nnew) = pack(source(1:nsrc), lsource_mask(1:nsrc)) + else + arr(nold + 1:nnew) = source(1:nsrc) + end if + + return + end subroutine base_util_append_arr_char_string + + + subroutine base_util_append_arr_DP(arr, source, nold, lsource_mask) + !! author: David A. Minton + !! + !! Append a single array of double precision type onto another. If the destination array is not allocated, or is not big enough, this will allocate space for it. + implicit none + ! Arguments + real(DP), dimension(:), allocatable, intent(inout) :: arr !! Destination array + real(DP), dimension(:), intent(in) :: source !! Array to append + integer(I4B), intent(in), optional :: nold !! Extent of original array. If passed, the source array will begin at arr(nold+1). Otherwise, the size of arr will be used. + logical, dimension(:), intent(in), optional :: lsource_mask !! Logical mask indicating which elements to append to + ! Internals + integer(I4B) :: nnew, nsrc, nend_orig + + if (present(lsource_mask)) then + nsrc = count(lsource_mask(:)) + else + nsrc = size(source) + end if + if (nsrc == 0) return + + if (.not.allocated(arr)) then + nend_orig = 0 + allocate(arr(nsrc)) + else + if (present(nold)) then + nend_orig = nold + else + nend_orig = size(arr) + end if + call util_resize(arr, nend_orig + nsrc) + end if + nnew = nend_orig + nsrc + + if (present(lsource_mask)) then + arr(nold + 1:nnew) = pack(source(1:nsrc), lsource_mask(1:nsrc)) + else + arr(nold + 1:nnew) = source(1:nsrc) + end if + + return + end subroutine base_util_append_arr_DP + + + subroutine base_util_append_arr_DPvec(arr, source, nold, lsource_mask) + !! author: David A. Minton + !! + !! Append a single array of double precision vector type of size (NDIM, n) onto another. If the destination array is not allocated, or is not big enough, this will allocate space for it. + implicit none + ! Arguments + real(DP), dimension(:,:), allocatable, intent(inout) :: arr !! Destination array + real(DP), dimension(:,:), intent(in) :: source !! Array to append + integer(I4B), intent(in), optional :: nold !! Extent of original array. If passed, the source array will begin at arr(nold+1). Otherwise, the size of arr will be used. + logical, dimension(:), intent(in), optional :: lsource_mask !! Logical mask indicating which elements to append to + ! Internals + integer(I4B) :: nnew, nsrc, nend_orig + + if (present(lsource_mask)) then + nsrc = count(lsource_mask(:)) + else + nsrc = size(source,dim=2) + end if + if (nsrc == 0) return + + if (.not.allocated(arr)) then + nend_orig = 0 + allocate(arr(NDIM,nsrc)) + else + if (present(nold)) then + nend_orig = nold + else + nend_orig = size(arr,dim=2) + end if + call util_resize(arr, nend_orig + nsrc) + end if + nnew = nend_orig + nsrc + + if (present(lsource_mask)) then + arr(1, nold + 1:nnew) = pack(source(1,1:nsrc), lsource_mask(1:nsrc)) + arr(2, nold + 1:nnew) = pack(source(2,1:nsrc), lsource_mask(1:nsrc)) + arr(3, nold + 1:nnew) = pack(source(3,1:nsrc), lsource_mask(1:nsrc)) + else + arr(:,nold + 1:nnew) = source(:,1:nsrc) + end if + + return + end subroutine base_util_append_arr_DPvec + + + subroutine base_util_append_arr_I4B(arr, source, nold, lsource_mask) + !! author: David A. Minton + !! + !! Append a single array of integer(I4B) onto another. If the destination array is not allocated, or is not big enough, this will allocate space for it. + implicit none + ! Arguments + integer(I4B), dimension(:), allocatable, intent(inout) :: arr !! Destination array + integer(I4B), dimension(:), intent(in) :: source !! Array to append + integer(I4B), intent(in), optional :: nold !! Extent of original array. If passed, the source array will begin at arr(nold+1). Otherwise, the size of arr will be used. + logical, dimension(:), intent(in), optional :: lsource_mask !! Logical mask indicating which elements to append to + ! Internals + integer(I4B) :: nnew, nsrc, nend_orig + + if (present(lsource_mask)) then + nsrc = count(lsource_mask(:)) + else + nsrc = size(source) + end if + if (nsrc == 0) return + + if (.not.allocated(arr)) then + nend_orig = 0 + allocate(arr(nsrc)) + else + if (present(nold)) then + nend_orig = nold + else + nend_orig = size(arr) + end if + call util_resize(arr, nend_orig + nsrc) + end if + nnew = nend_orig + nsrc + + if (present(lsource_mask)) then + arr(nold + 1:nnew) = pack(source(1:nsrc), lsource_mask(1:nsrc)) + else + arr(nold + 1:nnew) = source(1:nsrc) + end if + + return + end subroutine base_util_append_arr_I4B + + + subroutine base_util_append_arr_logical(arr, source, nold, lsource_mask) + !! author: David A. Minton + !! + !! Append a single array of logical type onto another. If the destination array is not allocated, or is not big enough, this will allocate space for it. + implicit none + ! Arguments + logical, dimension(:), allocatable, intent(inout) :: arr !! Destination array + logical, dimension(:), intent(in) :: source !! Array to append + integer(I4B), intent(in), optional :: nold !! Extent of original array. If passed, the source array will begin at arr(nold+1). Otherwise, the size of arr will be used. + logical, dimension(:), intent(in), optional :: lsource_mask !! Logical mask indicating which elements to append to + ! Internals + integer(I4B) :: nnew, nsrc, nend_orig + + if (present(lsource_mask)) then + nsrc = count(lsource_mask(:)) + else + nsrc = size(source) + end if + if (nsrc == 0) return + + if (.not.allocated(arr)) then + nend_orig = 0 + allocate(arr(nsrc)) + else + if (present(nold)) then + nend_orig = nold + else + nend_orig = size(arr) + end if + call util_resize(arr, nend_orig + nsrc) + end if + nnew = nend_orig + nsrc + + if (present(lsource_mask)) then + arr(nold + 1:nnew) = pack(source(1:nsrc), lsource_mask(1:nsrc)) + else + arr(nold + 1:nnew) = source(1:nsrc) + end if + + return + end subroutine base_util_append_arr_logical + + subroutine base_util_copy_store(self, source) !! author: David A. Minton !! @@ -303,6 +580,109 @@ subroutine base_util_exit(code) end subroutine base_util_exit + subroutine base_util_fill_arr_char_string(keeps, inserts, lfill_list) + !! author: David A. Minton + !! + !! Performs a fill operation on a single array of type character strings + !! This is the inverse of a spill operation + implicit none + ! Arguments + character(len=STRMAX), dimension(:), allocatable, intent(inout) :: keeps !! Array of values to keep + character(len=STRMAX), dimension(:), allocatable, intent(in) :: inserts !! Array of values to insert into keep + logical, dimension(:), intent(in) :: lfill_list !! Logical array of bodies to merge into the keeps + + if (.not.allocated(keeps) .or. .not.allocated(inserts)) return + + keeps(:) = unpack(keeps(:), .not.lfill_list(:), keeps(:)) + keeps(:) = unpack(inserts(:), lfill_list(:), keeps(:)) + + return + end subroutine base_util_fill_arr_char_string + + + subroutine base_util_fill_arr_DP(keeps, inserts, lfill_list) + !! author: David A. Minton + !! + !! Performs a fill operation on a single array of type DP + !! This is the inverse of a spill operation + implicit none + ! Arguments + real(DP), dimension(:), allocatable, intent(inout) :: keeps !! Array of values to keep + real(DP), dimension(:), allocatable, intent(in) :: inserts !! Array of values to insert into keep + logical, dimension(:), intent(in) :: lfill_list !! Logical array of bodies to merge into the keeps + + if (.not.allocated(keeps) .or. .not.allocated(inserts)) return + + keeps(:) = unpack(keeps(:), .not.lfill_list(:), keeps(:)) + keeps(:) = unpack(inserts(:), lfill_list(:), keeps(:)) + + return + end subroutine base_util_fill_arr_DP + + + subroutine base_util_fill_arr_DPvec(keeps, inserts, lfill_list) + !! author: David A. Minton + !! + !! Performs a fill operation on a single array of DP vectors with shape (NDIM, n) + !! This is the inverse of a spill operation + implicit none + ! Arguments + real(DP), dimension(:,:), allocatable, intent(inout) :: keeps !! Array of values to keep + real(DP), dimension(:,:), allocatable, intent(in) :: inserts !! Array of values to insert into keep + logical, dimension(:), intent(in) :: lfill_list !! Logical array of bodies to merge into the keeps + ! Internals + integer(I4B) :: i + + if (.not.allocated(keeps) .or. .not.allocated(inserts)) return + + do i = 1, NDIM + keeps(i,:) = unpack(keeps(i,:), .not.lfill_list(:), keeps(i,:)) + keeps(i,:) = unpack(inserts(i,:), lfill_list(:), keeps(i,:)) + end do + + return + end subroutine base_util_fill_arr_DPvec + + + subroutine base_util_fill_arr_I4B(keeps, inserts, lfill_list) + !! author: David A. Minton + !! + !! Performs a fill operation on a single array of type I4B + !! This is the inverse of a spill operation + implicit none + ! Arguments + integer(I4B), dimension(:), allocatable, intent(inout) :: keeps !! Array of values to keep + integer(I4B), dimension(:), allocatable, intent(in) :: inserts !! Array of values to insert into keep + logical, dimension(:), intent(in) :: lfill_list !! Logical array of bodies to merge into the keeps + + if (.not.allocated(keeps) .or. .not.allocated(inserts)) return + + keeps(:) = unpack(keeps(:), .not.lfill_list(:), keeps(:)) + keeps(:) = unpack(inserts(:), lfill_list(:), keeps(:)) + + return + end subroutine base_util_fill_arr_I4B + + + subroutine base_util_fill_arr_logical(keeps, inserts, lfill_list) + !! author: David A. Minton + !! + !! Performs a fill operation on a single array of logicals + !! This is the inverse of a spill operation + implicit none + ! Arguments + logical, dimension(:), allocatable, intent(inout) :: keeps !! Array of values to keep + logical, dimension(:), allocatable, intent(in) :: inserts !! Array of values to insert into keep + logical, dimension(:), intent(in) :: lfill_list !! Logical array of bodies to merge into the keeps + + if (.not.allocated(keeps) .or. .not.allocated(inserts)) return + + keeps(:) = unpack(keeps(:), .not.lfill_list(:), keeps(:)) + keeps(:) = unpack(inserts(:), lfill_list(:), keeps(:)) + + return + end subroutine base_util_fill_arr_logical + subroutine base_util_reset_storage(self) !! author: David A. Minton !! @@ -330,6 +710,237 @@ subroutine base_util_reset_storage(self) return end subroutine base_util_reset_storage + + subroutine base_util_resize_arr_char_string(arr, nnew) + !! author: David A. Minton + !! + !! Resizes an array component of type character string. nnew = 0 will deallocate. + implicit none + ! Arguments + character(len=STRMAX), dimension(:), allocatable, intent(inout) :: arr !! Array to resize + integer(I4B), intent(in) :: nnew !! New size + ! Internals + character(len=STRMAX), dimension(:), allocatable :: tmp !! Temporary storage array in case the input array is already allocated + integer(I4B) :: nold !! Old size + + if (nnew < 0) return + + if (nnew == 0) then + if (allocated(arr)) deallocate(arr) + return + end if + + if (allocated(arr)) then + nold = size(arr) + else + nold = 0 + end if + + if (nnew == nold) return + + allocate(tmp(nnew)) + if (nold > 0) then + if (nnew > nold) then + tmp(1:nold) = arr(1:nold) + tmp(nold+1:nnew) = "" + else + tmp(1:nnew) = arr(1:nnew) + end if + else + tmp(1:nnew) = "" + end if + call move_alloc(tmp, arr) + + return + end subroutine base_util_resize_arr_char_string + + + subroutine base_util_resize_arr_DP(arr, nnew) + !! author: David A. Minton + !! + !! Resizes an array component of double precision type. Passing nnew = 0 will deallocate. + implicit none + ! Arguments + real(DP), dimension(:), allocatable, intent(inout) :: arr !! Array to resize + integer(I4B), intent(in) :: nnew !! New size + ! Internals + real(DP), dimension(:), allocatable :: tmp !! Temporary storage array in case the input array is already allocated + integer(I4B) :: nold !! Old size + real(DP), parameter :: init_val = 0.0_DP + + if (nnew < 0) return + + if (nnew == 0) then + if (allocated(arr)) deallocate(arr) + return + end if + + if (allocated(arr)) then + nold = size(arr) + else + nold = 0 + end if + + if (nnew == nold) return + + allocate(tmp(nnew)) + if (nold > 0) then + if (nnew > nold) then + tmp(1:nold) = arr(1:nold) + tmp(nold+1:nnew) = init_val + else + tmp(1:nnew) = arr(1:nnew) + end if + else + tmp(1:nnew) = init_val + end if + call move_alloc(tmp, arr) + + return + end subroutine base_util_resize_arr_DP + + + subroutine base_util_resize_arr_DPvec(arr, nnew) + !! author: David A. Minton + !! + !! Resizes an array component of double precision vectors of size (NDIM, n). Passing nnew = 0 will deallocate. + implicit none + ! Arguments + real(DP), dimension(:,:), allocatable, intent(inout) :: arr !! Array to resize + integer(I4B), intent(in) :: nnew !! New size + ! Internals + real(DP), dimension(:,:), allocatable :: tmp !! Temporary storage array in case the input array is already allocated + integer(I4B) :: nold !! Old size + real(DP), dimension(NDIM), parameter :: init_val = 0.0_DP + integer(I4B) :: i + + if (nnew < 0) return + + if (nnew == 0) then + if (allocated(arr)) deallocate(arr) + return + end if + + if (allocated(arr)) then + nold = size(arr, dim=2) + else + nold = 0 + end if + + if (nnew == nold) return + + allocate(tmp(NDIM, nnew)) + if (nold > 0) then + if (nnew > nold) then + tmp(:,1:nold) = arr(:,1:nold) + do i = nold+1, nnew + tmp(:,i) = init_val(:) + end do + else + tmp(:,1:nnew) = arr(:,1:nnew) + end if + else + do i = 1, nnew + tmp(:, i) = init_val(:) + end do + end if + call move_alloc(tmp, arr) + + return + + return + end subroutine base_util_resize_arr_DPvec + + + subroutine base_util_resize_arr_I4B(arr, nnew) + !! author: David A. Minton + !! + !! Resizes an array component of integer type. Passing nnew = 0 will deallocate. + implicit none + ! Arguments + integer(I4B), dimension(:), allocatable, intent(inout) :: arr !! Array to resize + integer(I4B), intent(in) :: nnew !! New size + ! Internals + integer(I4B), dimension(:), allocatable :: tmp !! Temporary storage array in case the input array is already allocated + integer(I4B) :: nold !! Old size + integer(I4B), parameter :: init_val = -1 + + if (nnew < 0) return + + if (nnew == 0) then + if (allocated(arr)) deallocate(arr) + return + end if + + if (allocated(arr)) then + nold = size(arr) + else + nold = 0 + end if + + if (nnew == nold) return + + allocate(tmp(nnew)) + if (nold > 0) then + if (nnew > nold) then + tmp(1:nold) = arr(1:nold) + tmp(nold+1:nnew) = init_val + else + tmp(1:nnew) = arr(1:nnew) + end if + else + tmp(1:nnew) = init_val + end if + call move_alloc(tmp, arr) + + return + end subroutine base_util_resize_arr_I4B + + + subroutine base_util_resize_arr_logical(arr, nnew) + !! author: David A. Minton + !! + !! Resizes an array component of logical type. Passing nnew = 0 will deallocate. + implicit none + ! Arguments + logical, dimension(:), allocatable, intent(inout) :: arr !! Array to resize + integer(I4B), intent(in) :: nnew !! New size + ! Internals + logical, dimension(:), allocatable :: tmp !! Temporary storage array in case the input array is already allocated + integer(I4B) :: nold !! Old size + logical, parameter :: init_val = .false. + + if (nnew < 0) return + + if (nnew == 0) then + if (allocated(arr)) deallocate(arr) + return + end if + + if (allocated(arr)) then + nold = size(arr) + else + nold = 0 + end if + + if (nnew == nold) return + + allocate(tmp(nnew)) + if (nold > 0) then + if (nnew > nold) then + tmp(1:nold) = arr(1:nold) + tmp(nold+1:nnew) = init_val + else + tmp(1:nnew) = arr(1:nnew) + end if + else + tmp(1:nnew) = init_val + end if + call move_alloc(tmp, arr) + + return + end subroutine base_util_resize_arr_logical + subroutine base_util_resize_storage(self, nnew) !! author: David A. Minton @@ -405,35 +1016,1114 @@ subroutine base_util_snapshot_save(self, snapshot) return end subroutine base_util_snapshot_save - - subroutine base_final_storage(self) + + subroutine base_util_spill_arr_char_string(keeps, discards, lspill_list, ldestructive) !! author: David A. Minton !! - !! Finalizer for the storage object + !! Performs a spill operation on a single array of type character strings + !! This is the inverse of a spill operation implicit none ! Arguments - class(base_storage), intent(inout) :: self + character(len=STRMAX), dimension(:), allocatable, intent(inout) :: keeps !! Array of values to keep + character(len=STRMAX), dimension(:), allocatable, intent(inout) :: discards !! Array of discards + logical, dimension(:), intent(in) :: lspill_list !! Logical array of bodies to spill into the discards + logical, intent(in) :: ldestructive !! Logical flag indicating whether or not this operation should alter the keeps array or not + ! Internals + integer(I4B) :: nspill, nkeep, nlist + character(len=STRMAX), dimension(:), allocatable :: tmp !! Array of values to keep - call self%dealloc() - return - end subroutine base_final_storage + nkeep = count(.not.lspill_list(:)) + nspill = count(lspill_list(:)) + nlist = size(lspill_list(:)) + if (.not.allocated(keeps) .or. nspill == 0) return + if (.not.allocated(discards)) then + allocate(discards(nspill)) + else if (size(discards) /= nspill) then + deallocate(discards) + allocate(discards(nspill)) + end if + + discards(:) = pack(keeps(1:nlist), lspill_list(1:nlist)) + if (ldestructive) then + if (nkeep > 0) then + allocate(tmp(nkeep)) + tmp(:) = pack(keeps(1:nlist), .not. lspill_list(1:nlist)) + call move_alloc(tmp, keeps) + else + deallocate(keeps) + end if + end if - subroutine base_final_storage_frame(self) - !! author: David A. Minton - !! - !! Finalizer for the storage frame data type - implicit none - type(base_storage_frame) :: self - - if (allocated(self%item)) deallocate(self%item) - return - end subroutine base_final_storage_frame + end subroutine base_util_spill_arr_char_string + -#ifdef COARRAY - subroutine base_coclone_param(self) + subroutine base_util_spill_arr_DP(keeps, discards, lspill_list, ldestructive) + !! author: David A. Minton + !! + !! Performs a spill operation on a single array of type DP + !! This is the inverse of a spill operation + implicit none + ! Arguments + real(DP), dimension(:), allocatable, intent(inout) :: keeps !! Array of values to keep + real(DP), dimension(:), allocatable, intent(inout) :: discards !! Array of discards + logical, dimension(:), intent(in) :: lspill_list !! Logical array of bodies to spill into the discardss + logical, intent(in) :: ldestructive !! Logical flag indicating whether or not this operation should alter the keeps array or not + ! Internals + integer(I4B) :: nspill, nkeep, nlist + real(DP), dimension(:), allocatable :: tmp !! Array of values to keep + + nkeep = count(.not.lspill_list(:)) + nspill = count(lspill_list(:)) + nlist = size(lspill_list(:)) + + if (.not.allocated(keeps) .or. nspill == 0) return + if (.not.allocated(discards)) then + allocate(discards(nspill)) + else if (size(discards) /= nspill) then + deallocate(discards) + allocate(discards(nspill)) + end if + + discards(:) = pack(keeps(1:nlist), lspill_list(1:nlist)) + if (ldestructive) then + if (nkeep > 0) then + allocate(tmp(nkeep)) + tmp(:) = pack(keeps(1:nlist), .not. lspill_list(1:nlist)) + call move_alloc(tmp, keeps) + else + deallocate(keeps) + end if + end if + + return + end subroutine base_util_spill_arr_DP + + + subroutine base_util_spill_arr_DPvec(keeps, discards, lspill_list, ldestructive) + !! author: David A. Minton + !! + !! Performs a spill operation on a single array of DP vectors with shape (NDIM, n) + !! This is the inverse of a spill operation + implicit none + ! Arguments + real(DP), dimension(:,:), allocatable, intent(inout) :: keeps !! Array of values to keep + real(DP), dimension(:,:), allocatable, intent(inout) :: discards !! Array discards + logical, dimension(:), intent(in) :: lspill_list !! Logical array of bodies to spill into the discards + logical, intent(in) :: ldestructive !! Logical flag indicating whether or not this operation should alter the keeps array or not + ! Internals + integer(I4B) :: i, nspill, nkeep, nlist + real(DP), dimension(:,:), allocatable :: tmp !! Array of values to keep + + nkeep = count(.not.lspill_list(:)) + nspill = count(lspill_list(:)) + nlist = size(lspill_list(:)) + + if (.not.allocated(keeps) .or. nspill == 0) return + if (.not.allocated(discards)) then + allocate(discards(NDIM, nspill)) + else if (size(discards, dim=2) /= nspill) then + deallocate(discards) + allocate(discards(NDIM, nspill)) + end if + + do i = 1, NDIM + discards(i,:) = pack(keeps(i,1:nlist), lspill_list(1:nlist)) + end do + if (ldestructive) then + if (nkeep > 0) then + allocate(tmp(NDIM, nkeep)) + do i = 1, NDIM + tmp(i, :) = pack(keeps(i, 1:nlist), .not. lspill_list(1:nlist)) + end do + call move_alloc(tmp, keeps) + else + deallocate(keeps) + end if + end if + + return + end subroutine base_util_spill_arr_DPvec + + + subroutine base_util_spill_arr_I4B(keeps, discards, lspill_list, ldestructive) + !! author: David A. Minton + !! + !! Performs a spill operation on a single array of type I4B + !! This is the inverse of a spill operation + implicit none + ! Arguments + integer(I4B), dimension(:), allocatable, intent(inout) :: keeps !! Array of values to keep + integer(I4B), dimension(:), allocatable, intent(inout) :: discards !! Array of discards + logical, dimension(:), intent(in) :: lspill_list !! Logical array of bodies to spill into the discards + logical, intent(in) :: ldestructive !! Logical flag indicating whether or not this operation should alter the keeps array or not + ! Internals + integer(I4B) :: nspill, nkeep, nlist + integer(I4B), dimension(:), allocatable :: tmp !! Array of values to keep + + nkeep = count(.not.lspill_list(:)) + nspill = count(lspill_list(:)) + nlist = size(lspill_list(:)) + + if (.not.allocated(keeps) .or. nspill == 0) return + if (.not.allocated(discards)) then + allocate(discards(nspill)) + else if (size(discards) /= nspill) then + deallocate(discards) + allocate(discards(nspill)) + end if + + discards(:) = pack(keeps(1:nlist), lspill_list(1:nlist)) + if (ldestructive) then + if (nkeep > 0) then + allocate(tmp(nkeep)) + tmp(:) = pack(keeps(1:nlist), .not. lspill_list(1:nlist)) + call move_alloc(tmp, keeps) + else + deallocate(keeps) + end if + end if + + return + end subroutine base_util_spill_arr_I4B + + + subroutine base_util_spill_arr_I8B(keeps, discards, lspill_list, ldestructive) + !! author: David A. Minton + !! + !! Performs a spill operation on a single array of type I4B + !! This is the inverse of a spill operation + implicit none + ! Arguments + integer(I8B), dimension(:), allocatable, intent(inout) :: keeps !! Array of values to keep + integer(I8B), dimension(:), allocatable, intent(inout) :: discards !! Array of discards + logical, dimension(:), intent(in) :: lspill_list !! Logical array of bodies to spill into the discards + logical, intent(in) :: ldestructive !! Logical flag indicating whether or not this operation should alter the keeps array or not + ! Internals + integer(I4B) :: nspill, nkeep, nlist + integer(I8B), dimension(:), allocatable :: tmp !! Array of values to keep + + nkeep = count(.not.lspill_list(:)) + nspill = count(lspill_list(:)) + nlist = size(lspill_list(:)) + + if (.not.allocated(keeps) .or. nspill == 0) return + if (.not.allocated(discards)) then + allocate(discards(nspill)) + else if (size(discards) /= nspill) then + deallocate(discards) + allocate(discards(nspill)) + end if + + discards(:) = pack(keeps(1:nlist), lspill_list(1:nlist)) + if (ldestructive) then + if (nkeep > 0) then + allocate(tmp(nkeep)) + tmp(:) = pack(keeps(1:nlist), .not. lspill_list(1:nlist)) + call move_alloc(tmp, keeps) + else + deallocate(keeps) + end if + end if + + return + end subroutine base_util_spill_arr_I8B + + + subroutine base_util_spill_arr_logical(keeps, discards, lspill_list, ldestructive) + !! author: David A. Minton + !! + !! Performs a spill operation on a single array of logicals + !! This is the inverse of a spill operation + implicit none + ! Arguments + logical, dimension(:), allocatable, intent(inout) :: keeps !! Array of values to keep + logical, dimension(:), allocatable, intent(inout) :: discards !! Array of discards + logical, dimension(:), intent(in) :: lspill_list !! Logical array of bodies to spill into the discards + logical, intent(in) :: ldestructive !! Logical flag indicating whether or not this operation should alter the keeps array or no + ! Internals + integer(I4B) :: nspill, nkeep, nlist + logical, dimension(:), allocatable :: tmp !! Array of values to keep + + nkeep = count(.not.lspill_list(:)) + nspill = count(lspill_list(:)) + nlist = size(lspill_list(:)) + + if (.not.allocated(keeps) .or. nspill == 0) return + if (.not.allocated(discards)) then + allocate(discards(nspill)) + else if (size(discards) /= nspill) then + deallocate(discards) + allocate(discards(nspill)) + end if + + discards(:) = pack(keeps(1:nlist), lspill_list(1:nlist)) + if (ldestructive) then + if (nkeep > 0) then + allocate(tmp(nkeep)) + tmp(:) = pack(keeps(1:nlist), .not. lspill_list(1:nlist)) + call move_alloc(tmp, keeps) + else + deallocate(keeps) + end if + end if + + return + end subroutine base_util_spill_arr_logical + + + pure subroutine base_util_sort_dp(arr) + !! author: David A. Minton + !! + !! Sort input DP precision array in place into ascending numerical order using quicksort. + !! + implicit none + ! Arguments + real(DP), dimension(:), intent(inout) :: arr + + call base_util_sort_qsort_DP(arr) + + return + end subroutine base_util_sort_dp + + + pure subroutine base_util_sort_index_dp(arr, ind) + !! author: David A. Minton + !! + !! Sort input DP precision array by index in ascending numerical order using quick sort. + !! This algorithm works well for partially sorted arrays (which is usually the case here). + !! If ind is supplied already allocated, we assume it is an existing index array (e.g. a previously + !! sorted array). If it is not allocated, this subroutine swiftest_allocates it. + !! + implicit none + ! Arguments + real(DP), dimension(:), intent(in) :: arr + integer(I4B), dimension(:), allocatable, intent(inout) :: ind + ! Internals + integer(I4B) :: n, i + real(DP), dimension(:), allocatable :: tmparr + + n = size(arr) + if (.not.allocated(ind)) then + allocate(ind(n)) + ind = [(i, i=1, n)] + end if + allocate(tmparr, mold=arr) + tmparr(:) = arr(ind(:)) + call base_util_sort_qsort_DP(tmparr, ind) + + return + end subroutine base_util_sort_index_dp + + + recursive pure subroutine base_util_sort_qsort_DP(arr, ind) + !! author: David A. Minton + !! + !! Sort input DP precision array by index in ascending numerical order using quicksort sort. + !! + implicit none + ! Arguments + real(DP), dimension(:), intent(inout) :: arr + integer(I4B),dimension(:),intent(out), optional :: ind + !! Internals + integer :: iq + + if (size(arr) > 1) then + if (present(ind)) then + call base_util_sort_partition_DP(arr, iq, ind) + call base_util_sort_qsort_DP(arr(:iq-1),ind(:iq-1)) + call base_util_sort_qsort_DP(arr(iq:), ind(iq:)) + else + call base_util_sort_partition_DP(arr, iq) + call base_util_sort_qsort_DP(arr(:iq-1)) + call base_util_sort_qsort_DP(arr(iq:)) + end if + end if + + return + end subroutine base_util_sort_qsort_DP + + + pure subroutine base_util_sort_partition_DP(arr, marker, ind) + !! author: David A. Minton + !! + !! Partition function for quicksort on DP type + !! + implicit none + ! Arguments + real(DP), intent(inout), dimension(:) :: arr + integer(I4B), intent(inout), dimension(:), optional :: ind + integer(I4B), intent(out) :: marker + ! Internals + integer(I4B) :: i, j, itmp, narr, ipiv + real(DP) :: temp + real(DP) :: x ! pivot point + + narr = size(arr) + + ! Get center as pivot, as this is likely partially sorted + ipiv = narr / 2 + x = arr(ipiv) + i = 0 + j = narr + 1 + + do + j = j - 1 + do + if (arr(j) <= x) exit + j = j - 1 + end do + i = i + 1 + do + if (arr(i) >= x) exit + i = i + 1 + end do + if (i < j) then + ! exchange A(i) and A(j) + temp = arr(i) + arr(i) = arr(j) + arr(j) = temp + if (present(ind)) then + itmp = ind(i) + ind(i) = ind(j) + ind(j) = itmp + end if + else if (i == j) then + marker = i + 1 + return + else + marker = i + return + endif + end do + + return + end subroutine base_util_sort_partition_DP + + + pure subroutine base_util_sort_i4b(arr) + !! author: David A. Minton + !! + !! Sort input integer array in place into ascending numerical order using quick sort. + !! This algorithm works well for partially sorted arrays (which is usually the case here) + !! + implicit none + ! Arguments + integer(I4B), dimension(:), intent(inout) :: arr + + call base_util_sort_qsort_I4B(arr) + + return + end subroutine base_util_sort_i4b + + + pure subroutine base_util_sort_index_I4B(arr, ind) + !! author: David A. Minton + !! + !! Sort input integer array by index in ascending numerical order using quicksort. + !! If ind is supplied already allocated, we assume it is an existing index array (e.g. a previously + !! sorted array). If it is not allocated, this subroutine swiftest_allocates it. + !! + implicit none + ! Arguments + integer(I4B), dimension(:), intent(in) :: arr + integer(I4B), dimension(:), allocatable, intent(inout) :: ind + ! Internals + integer(I4B) :: n, i + integer(I4B), dimension(:), allocatable :: tmparr + + n = size(arr) + if (.not.allocated(ind)) then + allocate(ind(n)) + ind = [(i, i=1, n)] + end if + allocate(tmparr, mold=arr) + tmparr(:) = arr(ind(:)) + call base_util_sort_qsort_I4B(tmparr, ind) + + return + end subroutine base_util_sort_index_I4B + + + pure subroutine base_util_sort_index_I4B_I8Bind(arr, ind) + !! author: David A. Minton + !! + !! Sort input integer array by index in ascending numerical order using quicksort. + !! If ind is supplied already allocated, we assume it is an existing index array (e.g. a previously + !! sorted array). If it is not allocated, this subroutine swiftest_allocates it. + !! + implicit none + ! Arguments + integer(I4B), dimension(:), intent(in) :: arr + integer(I8B), dimension(:), allocatable, intent(inout) :: ind + ! Internals + integer(I8B) :: n, i + integer(I4B), dimension(:), allocatable :: tmparr + + n = size(arr) + if (.not.allocated(ind)) then + allocate(ind(n)) + ind = [(i, i=1_I8B, n)] + end if + allocate(tmparr, mold=arr) + tmparr(:) = arr(ind(:)) + call base_util_sort_qsort_I4B_I8Bind(tmparr, ind) + + return + end subroutine base_util_sort_index_I4B_I8Bind + + + pure subroutine base_util_sort_index_I8B_I8Bind(arr, ind) + !! author: David A. Minton + !! + !! Sort input integer array by index in ascending numerical order using quicksort. + !! If ind is supplied already allocated, we assume it is an existing index array (e.g. a previously + !! sorted array). If it is not allocated, this subroutine swiftest_allocates it. + !! + implicit none + ! Arguments + integer(I8B), dimension(:), intent(in) :: arr + integer(I8B), dimension(:), allocatable, intent(inout) :: ind + ! Internals + integer(I8B) :: n, i + integer(I8B), dimension(:), allocatable :: tmparr + + n = size(arr) + if (.not.allocated(ind)) then + allocate(ind(n)) + ind = [(i, i=1_I8B, n)] + end if + allocate(tmparr, mold=arr) + tmparr(:) = arr(ind(:)) + call base_util_sort_qsort_I8B_I8Bind(tmparr, ind) + + return + end subroutine base_util_sort_index_I8B_I8Bind + + + recursive pure subroutine base_util_sort_qsort_I4B(arr, ind) + !! author: David A. Minton + !! + !! Sort input I4B array by index in ascending numerical order using quicksort. + !! + implicit none + ! Arguments + integer(I4B), dimension(:), intent(inout) :: arr + integer(I4B), dimension(:), intent(out), optional :: ind + ! Internals + integer(I4B) :: iq + + if (size(arr) > 1) then + if (present(ind)) then + call base_util_sort_partition_I4B(arr, iq, ind) + call base_util_sort_qsort_I4B(arr(:iq-1),ind(:iq-1)) + call base_util_sort_qsort_I4B(arr(iq:), ind(iq:)) + else + call base_util_sort_partition_I4B(arr, iq) + call base_util_sort_qsort_I4B(arr(:iq-1)) + call base_util_sort_qsort_I4B(arr(iq:)) + end if + end if + + return + end subroutine base_util_sort_qsort_I4B + + + recursive pure subroutine base_util_sort_qsort_I4B_I8Bind(arr, ind) + !! author: David A. Minton + !! + !! Sort input I4B array by index in ascending numerical order using quicksort. + !! + implicit none + ! Arguments + integer(I4B), dimension(:), intent(inout) :: arr + integer(I8B), dimension(:), intent(out), optional :: ind + ! Internals + integer(I8B) :: iq + + if (size(arr) > 1_I8B) then + if (present(ind)) then + call base_util_sort_partition_I4B_I8Bind(arr, iq, ind) + call base_util_sort_qsort_I4B_I8Bind(arr(:iq-1_I8B),ind(:iq-1_I8B)) + call base_util_sort_qsort_I4B_I8Bind(arr(iq:), ind(iq:)) + else + call base_util_sort_partition_I4B_I8Bind(arr, iq) + call base_util_sort_qsort_I4B_I8Bind(arr(:iq-1_I8B)) + call base_util_sort_qsort_I4B_I8Bind(arr(iq:)) + end if + end if + + return + end subroutine base_util_sort_qsort_I4B_I8Bind + + + recursive pure subroutine base_util_sort_qsort_I8B_I8Bind(arr, ind) + !! author: David A. Minton + !! + !! Sort input I8B array by index in ascending numerical order using quicksort. + !! + implicit none + ! Arguments + integer(I8B), dimension(:), intent(inout) :: arr + integer(I8B), dimension(:), intent(out), optional :: ind + ! Internals + integer(I8B) :: iq + + if (size(arr) > 1_I8B) then + if (present(ind)) then + call base_util_sort_partition_I8B_I8Bind(arr, iq, ind) + call base_util_sort_qsort_I8B_I8Bind(arr(:iq-1_I8B),ind(:iq-1_I8B)) + call base_util_sort_qsort_I8B_I8Bind(arr(iq:), ind(iq:)) + else + call base_util_sort_partition_I8B_I8Bind(arr, iq) + call base_util_sort_qsort_I8B_I8Bind(arr(:iq-1_I8B)) + call base_util_sort_qsort_I8B_I8Bind(arr(iq:)) + end if + end if + + return + end subroutine base_util_sort_qsort_I8B_I8Bind + + + pure subroutine base_util_sort_partition_I4B(arr, marker, ind) + !! author: David A. Minton + !! + !! Partition function for quicksort on I4B type + !! + implicit none + ! Arguments + integer(I4B), intent(inout), dimension(:) :: arr + integer(I4B), intent(inout), dimension(:), optional :: ind + integer(I4B), intent(out) :: marker + ! Internals + integer(I4B) :: i, j, itmp, narr, ipiv + integer(I4B) :: temp + integer(I4B) :: x ! pivot point + + narr = size(arr) + + ! Get center as pivot, as this is likely partially sorted + ipiv = narr / 2 + x = arr(ipiv) + i = 0 + j = narr + 1 + + do + j = j - 1 + do + if (arr(j) <= x) exit + j = j - 1 + end do + i = i + 1 + do + if (arr(i) >= x) exit + i = i + 1 + end do + if (i < j) then + ! exchange A(i) and A(j) + temp = arr(i) + arr(i) = arr(j) + arr(j) = temp + if (present(ind)) then + itmp = ind(i) + ind(i) = ind(j) + ind(j) = itmp + end if + else if (i == j) then + marker = i + 1 + return + else + marker = i + return + endif + end do + + return + end subroutine base_util_sort_partition_I4B + + + pure subroutine base_util_sort_partition_I4B_I8Bind(arr, marker, ind) + !! author: David A. Minton + !! + !! Partition function for quicksort on I4B type + !! + implicit none + ! Arguments + integer(I4B), intent(inout), dimension(:) :: arr + integer(I8B), intent(inout), dimension(:), optional :: ind + integer(I8B), intent(out) :: marker + ! Internals + integer(I8B) :: i, j, itmp, narr, ipiv + integer(I4B) :: temp + integer(I8B) :: x ! pivot point + + narr = size(arr) + + ! Get center as pivot, as this is likely partially sorted + ipiv = narr / 2_I8B + x = arr(ipiv) + i = 0_I8B + j = narr + 1_I8B + + do + j = j - 1_I8B + do + if (arr(j) <= x) exit + j = j - 1_I8B + end do + i = i + 1_I8B + do + if (arr(i) >= x) exit + i = i + 1_I8B + end do + if (i < j) then + ! exchange A(i) and A(j) + temp = arr(i) + arr(i) = arr(j) + arr(j) = temp + if (present(ind)) then + itmp = ind(i) + ind(i) = ind(j) + ind(j) = itmp + end if + else if (i == j) then + marker = i + 1_I8B + return + else + marker = i + return + endif + end do + + return + end subroutine base_util_sort_partition_I4B_I8Bind + + + pure subroutine base_util_sort_partition_I8B_I8Bind(arr, marker, ind) + !! author: David A. Minton + !! + !! Partition function for quicksort on I8B type with I8B index + !! + implicit none + ! Arguments + integer(I8B), intent(inout), dimension(:) :: arr + integer(I8B), intent(inout), dimension(:), optional :: ind + integer(I8B), intent(out) :: marker + ! Internals + integer(I8B) :: i, j, itmp, narr, ipiv + integer(I8B) :: temp + integer(I8B) :: x ! pivot point + + narr = size(arr) + + ! Get center as pivot, as this is likely partially sorted + ipiv = narr / 2_I8B + x = arr(ipiv) + i = 0_I8B + j = narr + 1_I8B + + do + j = j - 1_I8B + do + if (arr(j) <= x) exit + j = j - 1_I8B + end do + i = i + 1_I8B + do + if (arr(i) >= x) exit + i = i + 1_I8B + end do + if (i < j) then + ! exchange A(i) and A(j) + temp = arr(i) + arr(i) = arr(j) + arr(j) = temp + if (present(ind)) then + itmp = ind(i) + ind(i) = ind(j) + ind(j) = itmp + end if + else if (i == j) then + marker = i + 1_I8B + return + else + marker = i + return + endif + end do + + return + end subroutine base_util_sort_partition_I8B_I8Bind + + + pure subroutine base_util_sort_sp(arr) + !! author: David A. Minton + !! + !! Sort input DP precision array in place into ascending numerical order using quicksort. + !! + implicit none + ! Arguments + real(SP), dimension(:), intent(inout) :: arr + + call base_util_sort_qsort_SP(arr) + + return + end subroutine base_util_sort_sp + + + pure subroutine base_util_sort_index_sp(arr, ind) + !! author: David A. Minton + !! + !! Sort input DP precision array by index in ascending numerical order using quicksort. + !! If ind is supplied already allocated, we assume it is an existing index array (e.g. a previously + !! sorted array). If it is not allocated, this subroutine swiftest_allocates it. + !! + implicit none + ! Arguments + real(SP), dimension(:), intent(in) :: arr + integer(I4B), dimension(:), allocatable, intent(inout) :: ind + ! Internals + integer(I4B) :: n, i + real(SP), dimension(:), allocatable :: tmparr + + n = size(arr) + if (.not.allocated(ind)) then + allocate(ind(n)) + ind = [(i, i=1, n)] + end if + allocate(tmparr, mold=arr) + tmparr(:) = arr(ind(:)) + call base_util_sort_qsort_SP(tmparr, ind) + + return + end subroutine base_util_sort_index_sp + + + recursive pure subroutine base_util_sort_qsort_SP(arr, ind) + !! author: David A. Minton + !! + !! Sort input DP precision array by index in ascending numerical order using quicksort. + !! + implicit none + ! Arguments + real(SP), dimension(:), intent(inout) :: arr + integer(I4B),dimension(:),intent(out), optional :: ind + !! Internals + integer :: iq + + if (size(arr) > 1) then + if (present(ind)) then + call base_util_sort_partition_SP(arr, iq, ind) + call base_util_sort_qsort_SP(arr(:iq-1),ind(:iq-1)) + call base_util_sort_qsort_SP(arr(iq:), ind(iq:)) + else + call base_util_sort_partition_SP(arr, iq) + call base_util_sort_qsort_SP(arr(:iq-1)) + call base_util_sort_qsort_SP(arr(iq:)) + end if + end if + + return + end subroutine base_util_sort_qsort_SP + + + pure subroutine base_util_sort_partition_SP(arr, marker, ind) + !! author: David A. Minton + !! + !! Partition function for quicksort on SP type + !! + implicit none + ! Arguments + real(SP), intent(inout), dimension(:) :: arr + integer(I4B), intent(inout), dimension(:), optional :: ind + integer(I4B), intent(out) :: marker + ! Internals + integer(I4B) :: i, j, itmp, narr, ipiv + real(SP) :: temp + real(SP) :: x ! pivot point + + narr = size(arr) + + ! Get center as pivot, as this is likely partially sorted + ipiv = narr / 2 + x = arr(ipiv) + i = 0 + j = narr + 1 + + do + j = j - 1 + do + if (arr(j) <= x) exit + j = j - 1 + end do + i = i + 1 + do + if (arr(i) >= x) exit + i = i + 1 + end do + if (i < j) then + ! exchange A(i) and A(j) + temp = arr(i) + arr(i) = arr(j) + arr(j) = temp + if (present(ind)) then + itmp = ind(i) + ind(i) = ind(j) + ind(j) = itmp + end if + else if (i == j) then + marker = i + 1 + return + else + marker = i + return + endif + end do + + return + end subroutine base_util_sort_partition_SP + + + pure subroutine base_util_sort_rearrange_arr_char_string(arr, ind, n) + !! author: David A. Minton + !! + !! Rearrange a single array of character string in-place from an index list. + implicit none + ! Arguments + character(len=STRMAX), dimension(:), allocatable, intent(inout) :: arr !! Destination array + integer(I4B), dimension(:), intent(in) :: ind !! Index to rearrange against + integer(I4B), intent(in) :: n !! Number of elements in arr and ind to rearrange + ! Internals + character(len=STRMAX), dimension(:), allocatable :: tmp !! Temporary copy of arry used during rearrange operation + + if (.not. allocated(arr) .or. n <= 0) return + allocate(tmp, mold=arr) + tmp(1:n) = arr(ind) + call move_alloc(tmp, arr) + + return + end subroutine base_util_sort_rearrange_arr_char_string + + + pure subroutine base_util_sort_rearrange_arr_DP(arr, ind, n) + !! author: David A. Minton + !! + !! Rearrange a single array of DP type in-place from an index list. + implicit none + ! Arguments + real(DP), dimension(:), allocatable, intent(inout) :: arr !! Destination array + integer(I4B), dimension(:), intent(in) :: ind !! Index to rearrange against + integer(I4B), intent(in) :: n !! Number of elements in arr and ind to rearrange + ! Internals + real(DP), dimension(:), allocatable :: tmp !! Temporary copy of array used during rearrange operation + + if (.not. allocated(arr) .or. n <= 0) return + allocate(tmp, mold=arr) + tmp(1:n) = arr(ind) + call move_alloc(tmp, arr) + + return + end subroutine base_util_sort_rearrange_arr_DP + + + pure subroutine base_util_sort_rearrange_arr_DPvec(arr, ind, n) + !! author: David A. Minton + !! + !! Rearrange a single array of (NDIM,n) DP-type vectors in-place from an index list. + implicit none + ! Arguments + real(DP), dimension(:,:), allocatable, intent(inout) :: arr !! Destination array + integer(I4B), dimension(:), intent(in) :: ind !! Index to rearrange against + integer(I4B), intent(in) :: n !! Number of elements in arr and ind to rearrange + ! Internals + real(DP), dimension(:,:), allocatable :: tmp !! Temporary copy of array used during rearrange operation + + if (.not. allocated(arr) .or. n <= 0) return + allocate(tmp, mold=arr) + tmp(:,1:n) = arr(:, ind) + call move_alloc(tmp, arr) + + return + end subroutine base_util_sort_rearrange_arr_DPvec + + + pure subroutine base_util_sort_rearrange_arr_I4B(arr, ind, n) + !! author: David A. Minton + !! + !! Rearrange a single array of integers in-place from an index list. + implicit none + ! Arguments + integer(I4B), dimension(:), allocatable, intent(inout) :: arr !! Destination array + integer(I4B), dimension(:), intent(in) :: ind !! Index to rearrange against + integer(I4B), intent(in) :: n !! Number of elements in arr and ind to rearrange + ! Internals + integer(I4B), dimension(:), allocatable :: tmp !! Temporary copy of array used during rearrange operation + + if (.not. allocated(arr) .or. n <= 0) return + allocate(tmp, mold=arr) + tmp(1:n) = arr(ind) + call move_alloc(tmp, arr) + + return + end subroutine base_util_sort_rearrange_arr_I4B + + pure subroutine base_util_sort_rearrange_arr_I4B_I8Bind(arr, ind, n) + !! author: David A. Minton + !! + !! Rearrange a single array of integers in-place from an index list. + implicit none + ! Arguments + integer(I4B), dimension(:), allocatable, intent(inout) :: arr !! Destination array + integer(I8B), dimension(:), intent(in) :: ind !! Index to rearrange against + integer(I8B), intent(in) :: n !! Number of elements in arr and ind to rearrange + ! Internals + integer(I4B), dimension(:), allocatable :: tmp !! Temporary copy of array used during rearrange operation + + if (.not. allocated(arr) .or. n <= 0_I8B) return + allocate(tmp, mold=arr) + tmp(1:n) = arr(ind) + call move_alloc(tmp, arr) + + return + end subroutine base_util_sort_rearrange_arr_I4B_I8Bind + + + pure subroutine base_util_sort_rearrange_arr_logical_I8Bind(arr, ind, n) + !! author: David A. Minton + !! + !! Rearrange a single array of logicals in-place from an index list. + implicit none + ! Arguments + logical, dimension(:), allocatable, intent(inout) :: arr !! Destination array + integer(I8B), dimension(:), intent(in) :: ind !! Index to rearrange against + integer(I8B), intent(in) :: n !! Number of elements in arr and ind to rearrange + ! Internals + logical, dimension(:), allocatable :: tmp !! Temporary copy of array used during rearrange operation + + if (.not. allocated(arr) .or. n <= 0) return + allocate(tmp, mold=arr) + tmp(1:n) = arr(ind) + call move_alloc(tmp, arr) + + return + end subroutine base_util_sort_rearrange_arr_logical_I8Bind + + + pure subroutine base_util_sort_rearrange_arr_logical(arr, ind, n) + !! author: David A. Minton + !! + !! Rearrange a single array of logicals in-place from an index list. + implicit none + ! Arguments + logical, dimension(:), allocatable, intent(inout) :: arr !! Destination array + integer(I4B), dimension(:), intent(in) :: ind !! Index to rearrange against + integer(I4B), intent(in) :: n !! Number of elements in arr and ind to rearrange + ! Internals + logical, dimension(:), allocatable :: tmp !! Temporary copy of array used during rearrange operation + + if (.not. allocated(arr) .or. n <= 0) return + allocate(tmp, mold=arr) + tmp(1:n) = arr(ind) + call move_alloc(tmp, arr) + + return + end subroutine base_util_sort_rearrange_arr_logical + + + subroutine base_util_unique_DP(input_array, output_array, index_map) + !! author: David A. Minton + !! + !! Takes an input unsorted integer array and returns a new array of sorted, unique values (DP version) + implicit none + ! Arguments + real(DP), dimension(:), intent(in) :: input_array !! Unsorted input array + real(DP), dimension(:), allocatable, intent(out) :: output_array !! Sorted array of unique values + integer(I4B), dimension(:), allocatable, intent(out) :: index_map !! An array of the same size as input_array that such that any for any index i, output_array(index_map(i)) = input_array(i) + ! Internals + real(DP), dimension(:), allocatable :: unique_array + integer(I4B) :: n + real(DP) :: lo, hi + + allocate(unique_array, mold=input_array) + allocate(index_map(size(input_array))) + lo = minval(input_array) - 1 + hi = maxval(input_array) + + n = 0 + do + n = n + 1 + lo = minval(input_array(:), mask=input_array(:) > lo) + unique_array(n) = lo + where(input_array(:) == lo) index_map(:) = n + if (lo >= hi) exit + enddo + allocate(output_array(n), source=unique_array(1:n)) + + return + end subroutine base_util_unique_DP + + + subroutine base_util_unique_I4B(input_array, output_array, index_map) + !! author: David A. Minton + !! + !! Takes an input unsorted integer array and returns a new array of sorted, unique values (I4B version) + implicit none + ! Arguments + integer(I4B), dimension(:), intent(in) :: input_array !! Unsorted input array + integer(I4B), dimension(:), allocatable, intent(out) :: output_array !! Sorted array of unique values + integer(I4B), dimension(:), allocatable, intent(out) :: index_map !! An array of the same size as input_array that such that any for any index i, output_array(index_map(i)) = input_array(i) + ! Internals + integer(I4B), dimension(:), allocatable :: unique_array + integer(I4B) :: n, lo, hi + + allocate(unique_array, mold=input_array) + allocate(index_map, mold=input_array) + lo = minval(input_array) - 1 + hi = maxval(input_array) + + n = 0 + do + n = n + 1 + lo = minval(input_array(:), mask=input_array(:) > lo) + unique_array(n) = lo + where(input_array(:) == lo) index_map(:) = n + if (lo >= hi) exit + enddo + allocate(output_array(n), source=unique_array(1:n)) + + return + end subroutine base_util_unique_I4B + + + subroutine base_final_storage(self) + !! author: David A. Minton + !! + !! Finalizer for the storage object + implicit none + ! Arguments + class(base_storage), intent(inout) :: self + + call self%dealloc() + return + end subroutine base_final_storage + + + subroutine base_final_storage_frame(self) + !! author: David A. Minton + !! + !! Finalizer for the storage frame data type + implicit none + type(base_storage_frame) :: self + + if (allocated(self%item)) deallocate(self%item) + + return + end subroutine base_final_storage_frame + +#ifdef COARRAY + subroutine base_coclone_param(self) !! author: David A. Minton !! !! Broadcasts the image 1 parameter to all other images in a parameter coarray @@ -524,7 +2214,6 @@ subroutine base_coclone_param(self) return end subroutine base_coclone_param - #endif diff --git a/src/coarray/coarray_clone.f90 b/src/coarray/coarray_clone.f90 new file mode 100644 index 000000000..4e9a34591 --- /dev/null +++ b/src/coarray/coarray_clone.f90 @@ -0,0 +1,494 @@ +!! Copyright 2023 - David Minton +!! This file is part of Swiftest. +!! Swiftest is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License +!! as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version. +!! Swiftest is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty +!! of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details. +!! You should have received a copy of the GNU General Public License along with Swiftest. +!! If not, see: https://www.gnu.org/licenses. + +submodule (coarray) s_coarray_clone + use swiftest +contains + + module subroutine coarray_component_clone_char(var,src_img) + !! author: David A. Minton + !! + !! Copies a component of a coarray derived type from the specified source image to the current local one. The default source image is 1 + !! character scalar version + implicit none + ! Arguments + character(len=*), intent(inout) :: var + integer(I4B), intent(in),optional :: src_img + ! Internals + character(len=STRMAX),allocatable :: tmp[:] + integer(I4B) :: img, si + + allocate(tmp[*]) + if (present(src_img)) then + si = src_img + else + si = 1 + end if + + sync all + if (this_image() == si) then + do img = 1, num_images() + tmp[img] = var + end do + sync images(*) + else + sync images(si) + var = tmp[si] + end if + + return + end subroutine coarray_component_clone_char + + + module subroutine coarray_component_clone_DP(var,src_img) + !! author: David A. Minton + !! + !! Copies a component of a coarray derived type from the specified source image to the current local one. The default source image is 1 + !! real(DP) scalar version + implicit none + ! Arguments + real(DP), intent(inout) :: var + integer(I4B), intent(in),optional :: src_img + ! Internals + real(DP),allocatable :: tmp[:] + integer(I4B) :: img, si + + allocate(tmp[*]) + + if (present(src_img)) then + si = src_img + else + si = 1 + end if + + sync all + if (this_image() == si) then + do img = 1, num_images() + tmp[img] = var + end do + sync images(*) + else + sync images(si) + var = tmp[si] + end if + + + return + end subroutine coarray_component_clone_DP + + + module subroutine coarray_component_clone_DP_arr1D(var,src_img) + !! author: David A. Minton + !! + !! Copies a component of a coarray derived type from the specified source image to the current local one. The default source image is 1 + !! real(DP) 1D allocatable array version + implicit none + ! Arguments + real(DP), dimension(:), allocatable, intent(inout) :: var + integer(I4B), intent(in),optional :: src_img + ! Internals + real(DP), dimension(:), codimension[:], allocatable :: tmp + integer(I4B) :: img, si + integer(I4B), allocatable :: n[:] + logical, allocatable :: isalloc[:] + + allocate(isalloc[*]) + allocate(n[*]) + + if (present(src_img)) then + si = src_img + else + si = 1 + end if + + isalloc = allocated(var) + if (isalloc) n = size(var) + sync all + if (.not. isalloc[si]) return + + allocate(tmp(n[si])[*]) + if (this_image() == si) then + do img = 1, num_images() + tmp(:)[img] = var + end do + sync images(*) + else + sync images(si) + if (allocated(var)) deallocate(var) + allocate(var, source=tmp) + end if + + return + end subroutine coarray_component_clone_DP_arr1D + + + module subroutine coarray_component_clone_DP_arr2D(var,src_img) + !! author: David A. Minton + !! + !! Copies a component of a coarray derived type from the specified source image to the current local one. The default source image is 1 + !! real(DP) 2D allocatable array version + implicit none + ! Arguments + real(DP), dimension(:,:), allocatable, intent(inout) :: var + integer(I4B), intent(in),optional :: src_img + ! Internals + real(DP), dimension(:,:), codimension[:], allocatable :: tmp + integer(I4B) :: img, si + integer(I4B), allocatable :: n1[:], n2[:] + logical, allocatable :: isalloc[:] + + allocate(n1[*]) + allocate(n2[*]) + allocate(isalloc[*]) + + if (present(src_img)) then + si = src_img + else + si = 1 + end if + + isalloc = allocated(var) + if (isalloc) then + n1 = size(var,dim=1) + n2 = size(var,dim=2) + end if + sync all + if (.not. isalloc[si]) return + + allocate(tmp(n1[si],n2[si])[*]) + if (this_image() == si) then + do img = 1, num_images() + tmp(:,:)[img] = var(:,:) + end do + sync images(*) + else + sync images(si) + if (allocated(var)) deallocate(var) + allocate(var, source=tmp) + end if + + return + end subroutine coarray_component_clone_DP_arr2D + + + module subroutine coarray_component_clone_DP_vec1D(var,src_img) + !! author: David A. Minton + !! + !! Copies a component of a coarray derived type from the specified source image to the current local one. The default source image is 1 + !! real(DP) 1D (NDIM) array version + implicit none + ! Arguments + real(DP), dimension(:), intent(inout) :: var + integer(I4B), intent(in),optional :: src_img + ! Internals + real(DP), dimension(:), codimension[:], allocatable :: tmp + integer(I4B) :: img, si + integer(I4B), allocatable :: n[:] + + allocate(n[*]) + if (present(src_img)) then + si = src_img + else + si = 1 + end if + + allocate(tmp(NDIM)[*]) + sync all + if (this_image() == si) then + do img = 1, num_images() + tmp(:)[img] = var(:) + end do + sync images(*) + else + sync images(si) + var(:) = tmp(:)[si] + end if + + deallocate(tmp) + + return + end subroutine coarray_component_clone_DP_vec1D + + + module subroutine coarray_component_clone_DP_vec2D(var,src_img) + !! author: David A. Minton + !! + !! Copies a component of a coarray derived type from the specified source image to the current local one. The default source image is 1 + !! real(DP) 1D allocatable array version + implicit none + ! Arguments + real(DP), dimension(:,:), allocatable, intent(inout) :: var + integer(I4B), intent(in),optional :: src_img + ! Internals + real(DP), dimension(:,:), codimension[:], allocatable :: tmp + integer(I4B) :: img, si + integer(I4B), allocatable :: n[:] + logical, allocatable :: isalloc[:] + + allocate(isalloc[*]) + allocate(n[*]) + + if (present(src_img)) then + si = src_img + else + si = 1 + end if + + isalloc = allocated(var) + if (isalloc) n = size(var,dim=2) + sync all + if (.not. isalloc[si]) return + + allocate(tmp(NDIM,n[si])[*]) + if (this_image() == si) then + do img = 1, num_images() + tmp(:,:)[img] = var(:,:) + end do + sync images(*) + else + sync images(si) + if (allocated(var)) deallocate(var) + allocate(var, source=tmp) + end if + + return + end subroutine coarray_component_clone_DP_vec2D + + + module subroutine coarray_component_clone_I4B(var,src_img) + !! author: David A. Minton + !! + !! Copies a component of a coarray derived type from the specified source image to the current local one. The default source image is 1 + !! integer(I4B) scalar version + implicit none + ! Arguments + integer(I4B), intent(inout) :: var + integer(I4B), intent(in),optional :: src_img + ! Internals + integer(I4B),allocatable :: tmp[:] + integer(I4B) :: img, si + + allocate(tmp[*]) + + if (present(src_img)) then + si = src_img + else + si = 1 + end if + + sync all + if (this_image() == si) then + do img = 1, num_images() + tmp[img] = var + end do + sync images(*) + else + sync images(si) + var = tmp[si] + end if + + return + end subroutine coarray_component_clone_I4B + + + module subroutine coarray_component_clone_I4B_arr1D(var,src_img) + !! author: David A. Minton + !! + !! Copies a component of a coarray derived type from the specified source image to the current local one. The default source image is 1 + !! integer(I4B) 1D allocatable array version + implicit none + ! Arguments + integer(I4B), dimension(:), allocatable, intent(inout) :: var + integer(I4B), intent(in),optional :: src_img + ! Internals + integer(I4B), dimension(:), codimension[:], allocatable :: tmp + integer(I4B) :: img, si + integer(I4B), allocatable :: n[:] + logical, allocatable :: isalloc[:] + + allocate(isalloc[*]) + allocate(n[*]) + if (present(src_img)) then + si = src_img + else + si = 1 + end if + + isalloc = allocated(var) + if (isalloc) n = size(var) + sync all + if (.not. isalloc[si]) return + + allocate(tmp(n[si])[*]) + if (this_image() == si) then + do img = 1, num_images() + tmp(:)[img] = var + end do + sync images(*) + else + sync images(si) + if (allocated(var)) deallocate(var) + allocate(var, source=tmp) + end if + + return + end subroutine coarray_component_clone_I4B_arr1D + + + module subroutine coarray_component_clone_I8B(var,src_img) + !! author: David A. Minton + !! + !! Copies a component of a coarray derived type from the specified source image to the current local one. The default source image is 1 + !! integer(I4B) scalar version + implicit none + ! Arguments + integer(I8B), intent(inout) :: var + integer(I4B), intent(in),optional :: src_img + ! Internals + integer(I8B),allocatable :: tmp[:] + integer(I4B) :: img, si + + allocate(tmp[*]) + if (present(src_img)) then + si = src_img + else + si = 1 + end if + + sync all + if (this_image() == si) then + do img = 1, num_images() + tmp[img] = var + end do + sync images(*) + else + sync images(si) + var = tmp[si] + end if + + return + end subroutine coarray_component_clone_I8B + + + module subroutine coarray_component_clone_lgt(var,src_img) + !! author: David A. Minton + !! + !! Copies a component of a coarray derived type from the specified source image to the current local one. The default source image is 1 + !! logical scalar version + implicit none + ! Arguments + logical, intent(inout) :: var + integer(I4B), intent(in),optional :: src_img + ! Internals + logical,allocatable :: tmp[:] + integer(I4B) :: img, si + + allocate(tmp[*]) + + if (present(src_img)) then + si = src_img + else + si = 1 + end if + + sync all + if (this_image() == si) then + do img = 1, num_images() + tmp[img] = var + end do + sync images(*) + else + sync images(si) + var = tmp[si] + end if + + + return + end subroutine coarray_component_clone_lgt + + + module subroutine coarray_component_clone_lgt_arr1D(var,src_img) + !! author: David A. Minton + !! + !! Copies a component of a coarray derived type from the specified source image to the current local one. The default source image is 1 + !! logical 1D allocatable array version + implicit none + ! Arguments + logical, dimension(:), allocatable, intent(inout) :: var + integer(I4B), intent(in),optional :: src_img + ! Internals + logical, dimension(:), codimension[:], allocatable :: tmp + integer(I4B) :: img, si + integer(I4B), allocatable :: n[:] + logical, allocatable :: isalloc[:] + + allocate(isalloc[*]) + allocate(n[*]) + if (present(src_img)) then + si = src_img + else + si = 1 + end if + + isalloc = allocated(var) + if (isalloc) n = size(var) + sync all + if (.not. isalloc[si]) return + + allocate(tmp(n[si])[*]) + if (this_image() == si) then + do img = 1, num_images() + tmp(:)[img] = var + end do + sync images(*) + else + sync images(si) + if (allocated(var)) deallocate(var) + allocate(var, source=tmp) + end if + + return + end subroutine coarray_component_clone_lgt_arr1D + + + module subroutine coarray_component_clone_QP(var,src_img) + !! author: David A. Minton + !! + !! Copies a component of a coarray derived type from the specified source image to the current local one. The default source image is 1 + !! real(DP) scalar version + implicit none + ! Arguments + real(QP), intent(inout) :: var + integer(I4B), intent(in),optional :: src_img + ! Internals + real(QP),allocatable :: tmp[:] + integer(I4B) :: img, si + + allocate(tmp[*]) + if (present(src_img)) then + si = src_img + else + si = 1 + end if + + !sync all + if (this_image() == si) then + do img = 1, num_images() + tmp[img] = var + end do + sync images(*) + else + sync images(si) + var = tmp[si] + end if + + return + end subroutine coarray_component_clone_QP + +end submodule s_coarray_clone \ No newline at end of file diff --git a/src/coarray/coarray_collect.f90 b/src/coarray/coarray_collect.f90 new file mode 100644 index 000000000..f9a65dad1 --- /dev/null +++ b/src/coarray/coarray_collect.f90 @@ -0,0 +1,296 @@ +!! Copyright 2023 - David Minton +!! This file is part of Swiftest. +!! Swiftest is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License +!! as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version. +!! Swiftest is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty +!! of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details. +!! You should have received a copy of the GNU General Public License along with Swiftest. +!! If not, see: https://www.gnu.org/licenses. + +submodule (coarray) s_coarray_collect + use swiftest +contains + + module subroutine coarray_component_collect_DP_arr1D(var,dest_img) + !! author: David A. Minton + !! + !! Collects components of a coarray derived type from all images and combines them into destination image component . The default destination image is 1 + !! real(DP) 1D allocatable array version + implicit none + ! Arguments + real(DP), dimension(:), allocatable, intent(inout) :: var + integer(I4B), intent(in),optional :: dest_img + ! Internals + real(DP), dimension(:), codimension[:], allocatable :: tmp + integer(I4B) :: i,img, ti, di, istart, iend, nmax + integer(I4B), allocatable :: n[:] + logical, allocatable :: isalloc[:] + + allocate(isalloc[*]) + allocate(n[*]) + + if (present(dest_img)) then + di = dest_img + else + di = 1 + end if + + isalloc = allocated(var) + if (isalloc) then + n = size(var) + else + n = 0 + end if + + sync all + + nmax = 0 + do img = 1, num_images() + if (n[img] > nmax) nmax = n[img] + end do + + allocate(tmp(nmax)[*]) + tmp(1:n) = var(1:n) + + !sync all + if (this_image() == di) then + do img = 1, num_images() + if (img /= di) then + call util_append(var, tmp(:)[img]) + n = n + n[img] + end if + end do + end if + + return + end subroutine coarray_component_collect_DP_arr1D + + + module subroutine coarray_component_collect_DP_arr2D(var,dest_img) + !! author: David A. Minton + !! + !! Collects components of a coarray derived type from all images and combines them into destination image component . The default destination image is 1 + !! real(DP) 2D allocatable array version + implicit none + ! Arguments + real(DP), dimension(:,:), allocatable, intent(inout) :: var + integer(I4B), intent(in),optional :: dest_img + ! Internals + real(DP), dimension(:,:), codimension[:], allocatable :: tmp + integer(I4B) :: i, img, ti, di, ntot, istart, iend, nmax + integer(I4B), allocatable :: n1[:], n2[:] + logical, allocatable :: isalloc[:] + + allocate(n1[*]) + allocate(n2[*]) + allocate(isalloc[*]) + + if (present(dest_img)) then + di = dest_img + else + di = 1 + end if + + isalloc = allocated(var) + if (isalloc) then + n1 = size(var,dim=1) + n2 = size(var,dim=2) + else + n1 = 0 + n2 = 0 + end if + sync all + + nmax = 0 + do img = 1, num_images() + if (n2[img] > nmax) nmax = n2[img] + end do + + allocate(tmp(NDIM,nmax)[*]) + tmp(:,1:n2) = var(:,1:n2) + + !sync all + if (this_image() == di) then + do img = 1, num_images() + if (img /= di) then + call util_append(var, tmp(:,:)[img]) + n2 = n2 + n2[img] + end if + end do + end if + + return + end subroutine coarray_component_collect_DP_arr2D + + + module subroutine coarray_component_collect_I4B(var,dest_img) + !! author: David A. Minton + !! + !! Sums this component of a coarray derived type from all images and places the value in the destination image component. The default destination image is 1 + !! integer(I4B) version + implicit none + ! Arguments + integer(I4B), intent(inout) :: var + integer(I4B), intent(in),optional :: dest_img + ! Internals + integer(I4B), allocatable :: tmp[:] + integer(I4B) :: img, di + + if (present(dest_img)) then + di = dest_img + else + di = 1 + end if + + allocate(tmp[*], source=var) + + !sync all + if (this_image() == di) then + var = 0 + do img = 1, num_images() + var = var + tmp[img] + end do + else + var = 0 + end if + + return + end subroutine coarray_component_collect_I4B + + + module subroutine coarray_component_collect_I8B(var,dest_img) + !! author: David A. Minton + !! + !! Sums this component of a coarray derived type from all images and places the value in the destination image component. The default destination image is 1 + !! integer(I8B) version + implicit none + ! Arguments + integer(I8B), intent(inout) :: var + integer(I4B), intent(in),optional :: dest_img + ! Internals + integer(I8B), allocatable :: tmp[:] + integer(I4B) :: img, di + + if (present(dest_img)) then + di = dest_img + else + di = 1 + end if + + allocate(tmp[*], source=var) + + !sync all + if (this_image() == di) then + var = 0 + do img = 1, num_images() + var = var + tmp[img] + end do + else + var = 0 + end if + + return + end subroutine coarray_component_collect_I8B + + + module subroutine coarray_component_collect_I4B_arr1D(var,dest_img) + !! author: David A. Minton + !! + !! Collects components of a coarray derived type from all images and combines them into destination image component . The default destination image is 1 + !! integer(I4B) 1D allocatable array version + implicit none + ! Arguments + integer(I4B), dimension(:), allocatable, intent(inout) :: var + integer(I4B), intent(in),optional :: dest_img + ! Internals + integer(I4B), dimension(:), codimension[:], allocatable :: tmp + integer(I4B) :: i,img, ti, di, istart, iend, nmax + integer(I4B), allocatable :: n[:] + logical, allocatable :: isalloc[:] + + allocate(isalloc[*]) + allocate(n[*]) + + if (present(dest_img)) then + di = dest_img + else + di = 1 + end if + + isalloc = allocated(var) + if (isalloc) then + n = size(var) + else + n = 0 + end if + sync all + nmax = 0 + do img = 1, num_images() + if (n[img] > nmax) nmax = n[img] + end do + + allocate(tmp(nmax)[*]) + tmp(1:n) = var(1:n) + + !sync all + if (this_image() == di) then + do img = 1, num_images() + if (img /= di) then + call util_append(var, tmp(:)[img]) + n = n + n[img] + end if + end do + end if + + return + end subroutine coarray_component_collect_I4B_arr1D + + + module subroutine coarray_component_collect_lgt_arr1D(var,dest_img) + !! author: David A. Minton + !! + !! Collects components of a coarray derived type from all images and combines them into destination image component . The default destination image is 1 + !! logical 1D allocatable array version + implicit none + ! Arguments + logical, dimension(:), allocatable, intent(inout) :: var + integer(I4B), intent(in),optional :: dest_img + ! Internals + logical, dimension(:), codimension[:], allocatable :: tmp + integer(I4B) :: i,img, ti, di, ntot, istart, iend + integer(I4B), allocatable :: n[:] + logical, allocatable :: isalloc[:] + + allocate(isalloc[*]) + allocate(n[*]) + + if (present(dest_img)) then + di = dest_img + else + di = 1 + end if + + isalloc = allocated(var) + if (isalloc) then + n = size(var) + else + n = 0 + end if + allocate(tmp[*],source=var) + sync all + + if (this_image() == di) then + do img = 1, num_images() + if (img /= di) then + call util_append(var, tmp(:)[img]) + n = n + n[img] + end if + end do + end if + + return + end subroutine coarray_component_collect_lgt_arr1D + + + +end submodule s_coarray_collect \ No newline at end of file diff --git a/src/coarray/coarray_module.f90 b/src/coarray/coarray_module.f90 index 597945761..b034f8786 100644 --- a/src/coarray/coarray_module.f90 +++ b/src/coarray/coarray_module.f90 @@ -8,871 +8,135 @@ !! If not, see: https://www.gnu.org/licenses. module coarray - !! author: David A. Minton - !! - !! Utilities that are used for coarray test particles - !! - use globals - implicit none - public + !! author: David A. Minton + !! + !! Utilities that are used for coarray test particles + !! + use globals + implicit none + public + + interface coclone + module subroutine coarray_component_clone_char(var,src_img) + implicit none + character(len=*), intent(inout) :: var + integer(I4B), intent(in),optional :: src_img + end subroutine coarray_component_clone_char + + module subroutine coarray_component_clone_DP(var,src_img) + implicit none + real(DP), intent(inout) :: var + integer(I4B), intent(in),optional :: src_img + end subroutine coarray_component_clone_DP + + module subroutine coarray_component_clone_DP_arr1D(var,src_img) + implicit none + real(DP), dimension(:), allocatable, intent(inout) :: var + integer(I4B), intent(in),optional :: src_img + end subroutine coarray_component_clone_DP_arr1D + + module subroutine coarray_component_clone_DP_arr2D(var,src_img) + implicit none + real(DP), dimension(:,:), allocatable, intent(inout) :: var + integer(I4B), intent(in),optional :: src_img + end subroutine coarray_component_clone_DP_arr2D + + module subroutine coarray_component_clone_I4B(var,src_img) + implicit none + integer(I4B), intent(inout) :: var + integer(I4B), intent(in),optional :: src_img + end subroutine coarray_component_clone_I4B + + module subroutine coarray_component_clone_I4B_arr1D(var,src_img) + implicit none + integer(I4B), dimension(:), allocatable, intent(inout) :: var + integer(I4B), intent(in),optional :: src_img + end subroutine coarray_component_clone_I4B_arr1D + + module subroutine coarray_component_clone_I8B(var,src_img) + implicit none + integer(I8B), intent(inout) :: var + integer(I4B), intent(in),optional :: src_img + end subroutine coarray_component_clone_I8B + + module subroutine coarray_component_clone_lgt(var,src_img) + implicit none + logical, intent(inout) :: var + integer(I4B), intent(in),optional :: src_img + end subroutine coarray_component_clone_lgt + + module subroutine coarray_component_clone_lgt_arr1D(var,src_img) + implicit none + logical, dimension(:), allocatable, intent(inout) :: var + integer(I4B), intent(in),optional :: src_img + end subroutine coarray_component_clone_lgt_arr1D + + module subroutine coarray_component_clone_QP(var,src_img) + implicit none + real(QP), intent(inout) :: var + integer(I4B), intent(in),optional :: src_img + end subroutine coarray_component_clone_QP + end interface + + + interface coclonevec + module subroutine coarray_component_clone_DP_vec1D(var,src_img) + implicit none + real(DP), dimension(:), intent(inout) :: var + integer(I4B), intent(in),optional :: src_img + end subroutine coarray_component_clone_DP_vec1D + + module subroutine coarray_component_clone_DP_vec2D(var,src_img) + implicit none + real(DP), dimension(:,:), allocatable, intent(inout) :: var + integer(I4B), intent(in),optional :: src_img + end subroutine coarray_component_clone_DP_vec2D + end interface coclonevec + + + interface cocollect + module subroutine coarray_component_collect_DP_arr1D(var,dest_img) + implicit none + real(DP), dimension(:), allocatable, intent(inout) :: var + integer(I4B), intent(in),optional :: dest_img + end subroutine coarray_component_collect_DP_arr1D + + + module subroutine coarray_component_collect_DP_arr2D(var,dest_img) + implicit none + real(DP), dimension(:,:), allocatable, intent(inout) :: var + integer(I4B), intent(in),optional :: dest_img + end subroutine coarray_component_collect_DP_arr2D + + + module subroutine coarray_component_collect_I4B(var,dest_img) + implicit none + integer(I4B), intent(inout) :: var + integer(I4B), intent(in),optional :: dest_img + end subroutine coarray_component_collect_I4B + + + module subroutine coarray_component_collect_I8B(var,dest_img) + implicit none + integer(I8B), intent(inout) :: var + integer(I4B), intent(in),optional :: dest_img + end subroutine coarray_component_collect_I8B + + + module subroutine coarray_component_collect_I4B_arr1D(var,dest_img) + implicit none + integer(I4B), dimension(:), allocatable, intent(inout) :: var + integer(I4B), intent(in),optional :: dest_img + end subroutine coarray_component_collect_I4B_arr1D + + + module subroutine coarray_component_collect_lgt_arr1D(var,dest_img) + implicit none + logical, dimension(:), allocatable, intent(inout) :: var + integer(I4B), intent(in),optional :: dest_img + end subroutine coarray_component_collect_lgt_arr1D + end interface - interface coclone - module procedure coarray_component_clone_char - module procedure coarray_component_clone_DP - module procedure coarray_component_clone_DP_arr1D - module procedure coarray_component_clone_DP_arr2D - module procedure coarray_component_clone_I4B - module procedure coarray_component_clone_I4B_arr1D - module procedure coarray_component_clone_I4B_arr2D - module procedure coarray_component_clone_I8B - module procedure coarray_component_clone_lgt - module procedure coarray_component_clone_lgt_arr1D - module procedure coarray_component_clone_QP - end interface - - interface cocollect - module procedure coarray_component_collect_DP_arr1D - module procedure coarray_component_collect_DP_arr2D - module procedure coarray_component_collect_I4B - module procedure coarray_component_collect_I4B_arr1D - module procedure coarray_component_collect_I4B_arr2D - module procedure coarray_component_collect_I8B - module procedure coarray_component_collect_lgt_arr1D - end interface - - contains - - subroutine coarray_component_clone_char(var,src_img) - !! author: David A. Minton - !! - !! Copies a component of a coarray derived type from the specified source image to the current local one. The default source image is 1 - !! character scalar version - implicit none - ! Arguments - character(len=*), intent(inout) :: var - integer(I4B), intent(in),optional :: src_img - ! Internals - character(len=STRMAX),save :: tmp[*] - integer(I4B) :: img, si - - if (present(src_img)) then - si = src_img - else - si = 1 - end if - - sync all - if (this_image() == si) then - do img = 1, num_images() - tmp[img] = var - end do - sync images(*) - else - sync images(si) - var = tmp[si] - end if - - return - end subroutine coarray_component_clone_char - - - subroutine coarray_component_clone_DP(var,src_img) - !! author: David A. Minton - !! - !! Copies a component of a coarray derived type from the specified source image to the current local one. The default source image is 1 - !! real(DP) scalar version - implicit none - ! Arguments - real(DP), intent(inout) :: var - integer(I4B), intent(in),optional :: src_img - ! Internals - real(DP),save :: tmp[*] - integer(I4B) :: img, si - - if (present(src_img)) then - si = src_img - else - si = 1 - end if - - sync all - if (this_image() == si) then - do img = 1, num_images() - tmp[img] = var - end do - sync images(*) - else - sync images(si) - var = tmp[si] - end if - - - return - end subroutine coarray_component_clone_DP - - - subroutine coarray_component_clone_DP_arr1D(var,src_img) - !! author: David A. Minton - !! - !! Copies a component of a coarray derived type from the specified source image to the current local one. The default source image is 1 - !! real(DP) 1D allocatable array version - implicit none - ! Arguments - real(DP), dimension(:), allocatable, intent(inout) :: var - integer(I4B), intent(in),optional :: src_img - ! Internals - real(DP), dimension(:), codimension[:], allocatable :: tmp - integer(I4B) :: img, si - integer(I4B), save :: n[*] - logical, save :: isalloc[*] - - if (present(src_img)) then - si = src_img - else - si = 1 - end if - - sync all - isalloc = allocated(var) - if (isalloc) n = size(var) - sync all - if (.not. isalloc[si]) return - - allocate(tmp(n[si])[*]) - if (this_image() == si) then - do img = 1, num_images() - tmp(:)[img] = var - end do - sync images(*) - else - sync images(si) - if (allocated(var)) deallocate(var) - allocate(var, source=tmp) - end if - - return - end subroutine coarray_component_clone_DP_arr1D - - - subroutine coarray_component_clone_DP_arr2D(var,src_img) - !! author: David A. Minton - !! - !! Copies a component of a coarray derived type from the specified source image to the current local one. The default source image is 1 - !! real(DP) 2D allocatable array version - implicit none - ! Arguments - real(DP), dimension(:,:), allocatable, intent(inout) :: var - integer(I4B), intent(in),optional :: src_img - ! Internals - real(DP), dimension(:,:), codimension[:], allocatable :: tmp - integer(I4B) :: img, si - integer(I4B), save :: n1[*], n2[*] - logical, save :: isalloc[*] - - if (present(src_img)) then - si = src_img - else - si = 1 - end if - - sync all - isalloc = allocated(var) - if (isalloc) then - n1 = size(var,dim=1) - n2 = size(var,dim=2) - end if - sync all - if (.not. isalloc[si]) return - - allocate(tmp(n1[si],n2[si])[*]) - if (this_image() == si) then - do img = 1, num_images() - tmp(:,:)[img] = var(:,:) - end do - sync images(*) - else - sync images(si) - if (allocated(var)) deallocate(var) - allocate(var, source=tmp) - end if - - return - end subroutine coarray_component_clone_DP_arr2D - - - subroutine coarray_component_clone_I4B(var,src_img) - !! author: David A. Minton - !! - !! Copies a component of a coarray derived type from the specified source image to the current local one. The default source image is 1 - !! integer(I4B) scalar version - implicit none - ! Arguments - integer(I4B), intent(inout) :: var - integer(I4B), intent(in),optional :: src_img - ! Internals - integer(I4B),save :: tmp[*] - integer(I4B) :: img, si - - if (present(src_img)) then - si = src_img - else - si = 1 - end if - - sync all - if (this_image() == si) then - do img = 1, num_images() - tmp[img] = var - end do - sync images(*) - else - sync images(si) - var = tmp[si] - end if - - return - end subroutine coarray_component_clone_I4B - - - subroutine coarray_component_clone_I4B_arr1D(var,src_img) - !! author: David A. Minton - !! - !! Copies a component of a coarray derived type from the specified source image to the current local one. The default source image is 1 - !! integer(I4B) 1D allocatable array version - implicit none - ! Arguments - integer(I4B), dimension(:), allocatable, intent(inout) :: var - integer(I4B), intent(in),optional :: src_img - ! Internals - integer(I4B), dimension(:), codimension[:], allocatable :: tmp - integer(I4B) :: img, si - integer(I4B), save :: n[*] - logical, save :: isalloc[*] - - if (present(src_img)) then - si = src_img - else - si = 1 - end if - - sync all - isalloc = allocated(var) - if (isalloc) n = size(var) - sync all - if (.not. isalloc[si]) return - - allocate(tmp(n[si])[*]) - if (this_image() == si) then - do img = 1, num_images() - tmp(:)[img] = var - end do - sync images(*) - else - sync images(si) - if (allocated(var)) deallocate(var) - allocate(var, source=tmp) - end if - - return - end subroutine coarray_component_clone_I4B_arr1D - - - subroutine coarray_component_clone_I4B_arr2D(var,src_img) - !! author: David A. Minton - !! - !! Copies a component of a coarray derived type from the specified source image to the current local one. The default source image is 1 - !! integer(I4B) 2D allocatable array version - implicit none - ! Arguments - integer(I4B), dimension(:,:), allocatable, intent(inout) :: var - integer(I4B), intent(in),optional :: src_img - ! Internals - integer(I4B), dimension(:,:), codimension[:], allocatable :: tmp - integer(I4B) :: img, si - integer(I4B), save :: n1[*], n2[*] - logical, save :: isalloc[*] - - if (present(src_img)) then - si = src_img - else - si = 1 - end if - - sync all - isalloc = allocated(var) - if (isalloc) then - n1 = size(var,dim=1) - n2 = size(var,dim=2) - end if - sync all - if (.not. isalloc[si]) return - - allocate(tmp(n1[si],n2[si])[*]) - if (this_image() == si) then - do img = 1, num_images() - tmp(:,:)[img] = var(:,:) - end do - sync images(*) - else - sync images(si) - if (allocated(var)) deallocate(var) - allocate(var, source=tmp) - end if - - return - end subroutine coarray_component_clone_I4B_arr2D - - - subroutine coarray_component_clone_I8B(var,src_img) - !! author: David A. Minton - !! - !! Copies a component of a coarray derived type from the specified source image to the current local one. The default source image is 1 - !! integer(I4B) scalar version - implicit none - ! Arguments - integer(I8B), intent(inout) :: var - integer(I4B), intent(in),optional :: src_img - ! Internals - integer(I8B),save :: tmp[*] - integer(I4B) :: img, si - - if (present(src_img)) then - si = src_img - else - si = 1 - end if - - sync all - if (this_image() == si) then - do img = 1, num_images() - tmp[img] = var - end do - sync images(*) - else - sync images(si) - var = tmp[si] - end if - - return - end subroutine coarray_component_clone_I8B - - - subroutine coarray_component_clone_lgt(var,src_img) - !! author: David A. Minton - !! - !! Copies a component of a coarray derived type from the specified source image to the current local one. The default source image is 1 - !! logical scalar version - implicit none - ! Arguments - logical, intent(inout) :: var - integer(I4B), intent(in),optional :: src_img - ! Internals - logical,save :: tmp[*] - integer(I4B) :: img, si - - if (present(src_img)) then - si = src_img - else - si = 1 - end if - - sync all - if (this_image() == si) then - do img = 1, num_images() - tmp[img] = var - end do - sync images(*) - else - sync images(si) - var = tmp[si] - end if - - - return - end subroutine coarray_component_clone_lgt - - - subroutine coarray_component_clone_lgt_arr1D(var,src_img) - !! author: David A. Minton - !! - !! Copies a component of a coarray derived type from the specified source image to the current local one. The default source image is 1 - !! logical 1D allocatable array version - implicit none - ! Arguments - logical, dimension(:), allocatable, intent(inout) :: var - integer(I4B), intent(in),optional :: src_img - ! Internals - logical, dimension(:), codimension[:], allocatable :: tmp - integer(I4B) :: img, si - integer(I4B), save :: n[*] - logical, save :: isalloc[*] - - if (present(src_img)) then - si = src_img - else - si = 1 - end if - - sync all - isalloc = allocated(var) - if (isalloc) n = size(var) - sync all - if (.not. isalloc[si]) return - - allocate(tmp(n[si])[*]) - if (this_image() == si) then - do img = 1, num_images() - tmp(:)[img] = var - end do - sync images(*) - else - sync images(si) - if (allocated(var)) deallocate(var) - allocate(var, source=tmp) - end if - - return - end subroutine coarray_component_clone_lgt_arr1D - - - subroutine coarray_component_clone_QP(var,src_img) - !! author: David A. Minton - !! - !! Copies a component of a coarray derived type from the specified source image to the current local one. The default source image is 1 - !! real(DP) scalar version - implicit none - ! Arguments - real(QP), intent(inout) :: var - integer(I4B), intent(in),optional :: src_img - ! Internals - real(QP),save :: tmp[*] - integer(I4B) :: img, si - - if (present(src_img)) then - si = src_img - else - si = 1 - end if - - sync all - if (this_image() == si) then - do img = 1, num_images() - tmp[img] = var - end do - sync images(*) - else - sync images(si) - var = tmp[si] - end if - - return - end subroutine coarray_component_clone_QP - - - subroutine coarray_component_collect_DP_arr1D(var,dest_img) - !! author: David A. Minton - !! - !! Collects components of a coarray derived type from all images and combines them into destination image component . The default destination image is 1 - !! real(DP) 1D allocatable array version - implicit none - ! Arguments - real(DP), dimension(:), allocatable, intent(inout) :: var - integer(I4B), intent(in),optional :: dest_img - ! Internals - real(DP), dimension(:), codimension[:], allocatable :: tmp - integer(I4B) :: img, ti, di, ntot, istart, iend - integer(I4B), save :: n[*] - logical, save :: isalloc[*] - - if (present(dest_img)) then - di = dest_img - else - di = 1 - end if - - sync all - isalloc = allocated(var) - if (isalloc) then - n = size(var) - else - n = 0 - end if - sync all - ntot = 0 - do img = 1, num_images() - ntot = ntot + n[img] - end do - - allocate(tmp(ntot)[*]) - - ti = this_image() - - istart = 1 - iend = n - do img = 1, this_image() - 1 - istart = istart + n[img] - iend = iend + n[img] - end do - - if (isalloc) then - tmp(istart:iend) = var(:) - deallocate(var) - end if - - sync all - if (this_image() == di) then - allocate(var(ntot)) - istart = 1 - iend = n - do img = 1, num_images() - var(istart:iend) = tmp[img](istart:iend) - istart = istart + n[img] - iend = iend + n[img] - end do - end if - - return - end subroutine coarray_component_collect_DP_arr1D - - - subroutine coarray_component_collect_DP_arr2D(var,dest_img) - !! author: David A. Minton - !! - !! Collects components of a coarray derived type from all images and combines them into destination image component . The default destination image is 1 - !! real(DP) 2D allocatable array version - implicit none - ! Arguments - real(DP), dimension(:,:), allocatable, intent(inout) :: var - integer(I4B), intent(in),optional :: dest_img - ! Internals - integer(I4B), dimension(:,:), codimension[:], allocatable :: tmp - integer(I4B) :: img, ti, di, ntot, istart, iend - integer(I4B), save :: n1[*], n2[*] - logical, save :: isalloc[*] - - if (present(dest_img)) then - di = dest_img - else - di = 1 - end if - - sync all - isalloc = allocated(var) - if (isalloc) then - n1 = size(var,dim=1) - n2 = size(var,dim=2) - else - n1 = 0 - n2 = 0 - end if - sync all - ntot = 0 - do img = 1, num_images() - ntot = ntot + n2[img] - end do - - allocate(tmp(n1,ntot)[*]) - - ti = this_image() - - istart = 1 - iend = n2 - do img = 1, this_image() - 1 - istart = istart + n2[img] - iend = iend + n2[img] - end do - - if (isalloc) then - tmp(:,istart:iend) = var(:,:) - deallocate(var) - end if - - sync all - if (this_image() == di) then - allocate(var(n1,ntot)) - - istart = 1 - iend = n2 - do img = 1, num_images() - var(:,istart:iend) = tmp[img](:,istart:iend) - istart = istart + n2[img] - iend = iend + n2[img] - end do - end if - - return - end subroutine coarray_component_collect_DP_arr2D - - - subroutine coarray_component_collect_I4B(var,dest_img) - !! author: David A. Minton - !! - !! Sums this component of a coarray derived type from all images and places the value in the destination image component. The default destination image is 1 - !! integer(I4B) version - implicit none - ! Arguments - integer(I4B), intent(inout) :: var - integer(I4B), intent(in),optional :: dest_img - ! Internals - integer(I4B), allocatable :: tmp[:] - integer(I4B) :: img, di - - if (present(dest_img)) then - di = dest_img - else - di = 1 - end if - - allocate(tmp[*], source=var) - sync all - - if (this_image() == di) then - var = 0 - do img = 1, num_images() - var = var + tmp[img] - end do - else - var = 0 - end if - - return - end subroutine coarray_component_collect_I4B - - - subroutine coarray_component_collect_I8B(var,dest_img) - !! author: David A. Minton - !! - !! Sums this component of a coarray derived type from all images and places the value in the destination image component. The default destination image is 1 - !! integer(I8B) version - implicit none - ! Arguments - integer(I8B), intent(inout) :: var - integer(I4B), intent(in),optional :: dest_img - ! Internals - integer(I8B), allocatable :: tmp[:] - integer(I4B) :: img, di - - if (present(dest_img)) then - di = dest_img - else - di = 1 - end if - - allocate(tmp[*], source=var) - sync all - - if (this_image() == di) then - var = 0 - do img = 1, num_images() - var = var + tmp[img] - end do - else - var = 0 - end if - - return - end subroutine coarray_component_collect_I8B - - - subroutine coarray_component_collect_I4B_arr1D(var,dest_img) - !! author: David A. Minton - !! - !! Collects components of a coarray derived type from all images and combines them into destination image component . The default destination image is 1 - !! integer(I4B) 1D allocatable array version - implicit none - ! Arguments - integer(I4B), dimension(:), allocatable, intent(inout) :: var - integer(I4B), intent(in),optional :: dest_img - ! Internals - integer(I4B), dimension(:), codimension[:], allocatable :: tmp - integer(I4B) :: img, ti, di, ntot, istart, iend - integer(I4B), save :: n[*] - logical, save :: isalloc[*] - - if (present(dest_img)) then - di = dest_img - else - di = 1 - end if - - sync all - isalloc = allocated(var) - if (isalloc) then - n = size(var) - else - n = 0 - end if - sync all - ntot = 0 - do img = 1, num_images() - ntot = ntot + n[img] - end do - - allocate(tmp(ntot)[*]) - - ti = this_image() - - istart = 1 - iend = n - do img = 1, this_image() - 1 - istart = istart + n[img] - iend = iend + n[img] - end do - - if (isalloc) then - tmp(istart:iend) = var(:) - deallocate(var) - end if - - sync all - if (this_image() == di) then - allocate(var(ntot)) - istart = 1 - iend = n - do img = 1, num_images() - var(istart:iend) = tmp[img](istart:iend) - istart = istart + n[img] - iend = iend + n[img] - end do - end if - - return - end subroutine coarray_component_collect_I4B_arr1D - - - subroutine coarray_component_collect_I4B_arr2D(var,dest_img) - !! author: David A. Minton - !! - !! Collects components of a coarray derived type from all images and combines them into destination image component . The default destination image is 1 - !! integer(I4B) 2D allocatable array version - implicit none - ! Arguments - integer(I4B), dimension(:,:), allocatable, intent(inout) :: var - integer(I4B), intent(in),optional :: dest_img - ! Internals - integer(I4B), dimension(:,:), codimension[:], allocatable :: tmp - integer(I4B) :: img, ti, di, ntot, istart, iend - integer(I4B), save :: n1[*], n2[*] - logical, save :: isalloc[*] - - if (present(dest_img)) then - di = dest_img - else - di = 1 - end if - - sync all - isalloc = allocated(var) - if (isalloc) then - n1 = size(var,dim=1) - n2 = size(var,dim=2) - else - n1 = 0 - n2 = 0 - end if - sync all - ntot = 0 - do img = 1, num_images() - ntot = ntot + n2[img] - end do - - allocate(tmp(n1,ntot)[*]) - - ti = this_image() - - istart = 1 - iend = n2 - do img = 1, this_image() - 1 - istart = istart + n2[img] - iend = iend + n2[img] - end do - - if (isalloc) then - tmp(:,istart:iend) = var(:,:) - deallocate(var) - end if - - sync all - if (this_image() == di) then - allocate(var(n1,ntot)) - - istart = 1 - iend = n2 - do img = 1, num_images() - var(:,istart:iend) = tmp[img](:,istart:iend) - istart = istart + n2[img] - iend = iend + n2[img] - end do - end if - - return - end subroutine coarray_component_collect_I4B_arr2D - - - subroutine coarray_component_collect_lgt_arr1D(var,dest_img) - !! author: David A. Minton - !! - !! Collects components of a coarray derived type from all images and combines them into destination image component . The default destination image is 1 - !! logical 1D allocatable array version - implicit none - ! Arguments - logical, dimension(:), allocatable, intent(inout) :: var - integer(I4B), intent(in),optional :: dest_img - ! Internals - logical, dimension(:), codimension[:], allocatable :: tmp - integer(I4B) :: img, ti, di, ntot, istart, iend - integer(I4B), save :: n[*] - logical, save :: isalloc[*] - - if (present(dest_img)) then - di = dest_img - else - di = 1 - end if - - sync all - isalloc = allocated(var) - if (isalloc) then - n = size(var) - else - n = 0 - end if - sync all - ntot = 0 - do img = 1, num_images() - ntot = ntot + n[img] - end do - - allocate(tmp(ntot)[*]) - - ti = this_image() - - istart = 1 - iend = n - do img = 1, this_image() - 1 - istart = istart + n[img] - iend = iend + n[img] - end do - - if (isalloc) then - tmp(istart:iend) = var(:) - deallocate(var) - end if - - sync all - if (this_image() == di) then - allocate(var(ntot)) - istart = 1 - iend = n - do img = 1, num_images() - var(istart:iend) = tmp[img](istart:iend) - istart = istart + n[img] - iend = iend + n[img] - end do - end if - - return - end subroutine coarray_component_collect_lgt_arr1D end module coarray \ No newline at end of file diff --git a/src/collision/collision_util.f90 b/src/collision/collision_util.f90 index 9eeb1746f..95152c0aa 100644 --- a/src/collision/collision_util.f90 +++ b/src/collision/collision_util.f90 @@ -227,7 +227,7 @@ module subroutine collision_util_index_map(self) call self%get_index_values(idvals, tvals) ! Consolidate ids to only unique values - call swiftest_util_unique(idvals,self%idvals,self%idmap) + call util_unique(idvals,self%idvals,self%idmap) self%nid = size(self%idvals) ! Don't consolidate time values (multiple collisions can happen in a single time step) diff --git a/src/encounter/encounter_check.f90 b/src/encounter/encounter_check.f90 index c2013ce6b..06d0a76d4 100644 --- a/src/encounter/encounter_check.f90 +++ b/src/encounter/encounter_check.f90 @@ -103,10 +103,10 @@ module subroutine encounter_check_all_plplm(param, nplm, nplt, rplm, vplm, rplt, call move_alloc(ltmp, lvdotr) nenc = nenc + plmplt_nenc - call swiftest_util_sort(index1, ind) - call swiftest_util_sort_rearrange(index1, ind, nenc) - call swiftest_util_sort_rearrange(index2, ind, nenc) - call swiftest_util_sort_rearrange(lvdotr, ind, nenc) + call util_sort(index1, ind) + call util_sort_rearrange(index1, ind, nenc) + call util_sort_rearrange(index2, ind, nenc) + call util_sort_rearrange(lvdotr, ind, nenc) end if @@ -677,10 +677,10 @@ subroutine encounter_check_remove_duplicates(n, nenc, index1, index2, lvdotr) return end if - call swiftest_util_sort(index1, ind) - call swiftest_util_sort_rearrange(index1, ind, nenc) - call swiftest_util_sort_rearrange(index2, ind, nenc) - call swiftest_util_sort_rearrange(lvdotr, ind, nenc) + call util_sort(index1, ind) + call util_sort_rearrange(index1, ind, nenc) + call util_sort_rearrange(index2, ind, nenc) + call util_sort_rearrange(lvdotr, ind, nenc) ! Get the bounds on the bodies in the first index ibeg(:) = n @@ -706,7 +706,7 @@ subroutine encounter_check_remove_duplicates(n, nenc, index1, index2, lvdotr) khi = iend(i) nenci = khi - klo + 1_I8B if (allocated(ind)) deallocate(ind) - call swiftest_util_sort(index2(klo:khi), ind) + call util_sort(index2(klo:khi), ind) index2(klo:khi) = itmp(klo - 1_I8B + ind(:)) do j = klo + 1_I8B, khi if (index2(j) == index2(j - 1_I8B)) lencounter(j) = .false. @@ -746,7 +746,7 @@ pure module subroutine encounter_check_sort_aabb_1D(self, n, extent_arr) ! Internals integer(I8B) :: i, k - call swiftest_util_sort(extent_arr, self%ind) + call util_sort(extent_arr, self%ind) do concurrent(k = 1_I8B:2_I8B * n) i = self%ind(k) diff --git a/src/encounter/encounter_module.f90 b/src/encounter/encounter_module.f90 index f3c24c763..bcae2bd15 100644 --- a/src/encounter/encounter_module.f90 +++ b/src/encounter/encounter_module.f90 @@ -259,9 +259,9 @@ end subroutine encounter_util_setup_list module subroutine encounter_util_append_list(self, source, lsource_mask) implicit none - class(encounter_list), intent(inout) :: self !! Swiftest encounter list object - class(encounter_list), intent(in) :: source !! Source object to append - logical, dimension(:), intent(in) :: lsource_mask !! Logical mask indicating which elements to append to + class(encounter_list), intent(inout) :: self !! Swiftest encounter list object + class(encounter_list), intent(in) :: source !! Source object to append + logical, dimension(:), intent(in) :: lsource_mask !! Logical mask indicating which elements to append to end subroutine encounter_util_append_list module subroutine encounter_util_copy_list(self, source) diff --git a/src/encounter/encounter_util.f90 b/src/encounter/encounter_util.f90 index 01848d571..f935afe5a 100644 --- a/src/encounter/encounter_util.f90 +++ b/src/encounter/encounter_util.f90 @@ -18,28 +18,27 @@ module subroutine encounter_util_append_list(self, source, lsource_mask) !! This method will automatically resize the destination body if it is too small implicit none ! Arguments - class(encounter_list), intent(inout) :: self !! Swiftest encounter list object - class(encounter_list), intent(in) :: source !! Source object to append - logical, dimension(:), intent(in) :: lsource_mask !! Logical mask indicating which elements to append to + class(encounter_list), intent(inout) :: self !! Swiftest encounter list object + class(encounter_list), intent(in) :: source !! Source object to append + logical, dimension(:), intent(in) :: lsource_mask !! Logical mask indicating which elements to append to ! Internals integer(I4B) :: nold, nsrc nold = int(self%nenc, kind=I4B) - nsrc = int(source%nenc, kind=I4B) - call swiftest_util_append(self%tcollision, source%tcollision, nold, nsrc, lsource_mask) - call swiftest_util_append(self%lclosest, source%lclosest, nold, nsrc, lsource_mask) - call swiftest_util_append(self%lvdotr, source%lvdotr, nold, nsrc, lsource_mask) - call swiftest_util_append(self%status, source%status, nold, nsrc, lsource_mask) - call swiftest_util_append(self%index1, source%index1, nold, nsrc, lsource_mask) - call swiftest_util_append(self%index2, source%index2, nold, nsrc, lsource_mask) - call swiftest_util_append(self%id1, source%id1, nold, nsrc, lsource_mask) - call swiftest_util_append(self%id2, source%id2, nold, nsrc, lsource_mask) - call swiftest_util_append(self%r1, source%r1, nold, nsrc, lsource_mask) - call swiftest_util_append(self%r2, source%r2, nold, nsrc, lsource_mask) - call swiftest_util_append(self%v1, source%v1, nold, nsrc, lsource_mask) - call swiftest_util_append(self%v2, source%v2, nold, nsrc, lsource_mask) - call swiftest_util_append(self%level, source%level, nold, nsrc, lsource_mask) - self%nenc = nold + count(lsource_mask(1:nsrc)) + call util_append(self%tcollision, source%tcollision, nold, lsource_mask) + call util_append(self%lclosest, source%lclosest, nold, lsource_mask) + call util_append(self%lvdotr, source%lvdotr, nold, lsource_mask) + call util_append(self%status, source%status, nold, lsource_mask) + call util_append(self%index1, source%index1, nold, lsource_mask) + call util_append(self%index2, source%index2, nold, lsource_mask) + call util_append(self%id1, source%id1, nold, lsource_mask) + call util_append(self%id2, source%id2, nold, lsource_mask) + call util_append(self%r1, source%r1, nold, lsource_mask) + call util_append(self%r2, source%r2, nold, lsource_mask) + call util_append(self%v1, source%v1, nold, lsource_mask) + call util_append(self%v2, source%v2, nold, lsource_mask) + call util_append(self%level, source%level, nold, lsource_mask) + self%nenc = nold + count(lsource_mask(:)) return end subroutine encounter_util_append_list @@ -283,11 +282,11 @@ module subroutine encounter_util_index_map(self) call encounter_util_get_vals_storage(self, idvals, tvals) ! Consolidate ids to only unique values - call swiftest_util_unique(idvals,self%idvals,self%idmap) + call util_unique(idvals,self%idvals,self%idmap) self%nid = size(self%idvals) ! Consolidate time values to only unique values - call swiftest_util_unique(tvals,self%tvals,self%tmap) + call util_unique(tvals,self%tvals,self%tmap) self%nt = size(self%tvals) return @@ -431,19 +430,19 @@ module subroutine encounter_util_spill_list(self, discards, lspill_list, ldestru integer(I8B) :: nenc_old associate(keeps => self) - call swiftest_util_spill(keeps%tcollision, discards%tcollision, lspill_list, ldestructive) - call swiftest_util_spill(keeps%lvdotr, discards%lvdotr, lspill_list, ldestructive) - call swiftest_util_spill(keeps%lclosest, discards%lclosest, lspill_list, ldestructive) - call swiftest_util_spill(keeps%status, discards%status, lspill_list, ldestructive) - call swiftest_util_spill(keeps%index1, discards%index1, lspill_list, ldestructive) - call swiftest_util_spill(keeps%index2, discards%index2, lspill_list, ldestructive) - call swiftest_util_spill(keeps%id1, discards%id1, lspill_list, ldestructive) - call swiftest_util_spill(keeps%id2, discards%id2, lspill_list, ldestructive) - call swiftest_util_spill(keeps%r1, discards%r1, lspill_list, ldestructive) - call swiftest_util_spill(keeps%r2, discards%r2, lspill_list, ldestructive) - call swiftest_util_spill(keeps%v1, discards%v1, lspill_list, ldestructive) - call swiftest_util_spill(keeps%v2, discards%v2, lspill_list, ldestructive) - call swiftest_util_spill(keeps%level, discards%level, lspill_list, ldestructive) + call util_spill(keeps%tcollision, discards%tcollision, lspill_list, ldestructive) + call util_spill(keeps%lvdotr, discards%lvdotr, lspill_list, ldestructive) + call util_spill(keeps%lclosest, discards%lclosest, lspill_list, ldestructive) + call util_spill(keeps%status, discards%status, lspill_list, ldestructive) + call util_spill(keeps%index1, discards%index1, lspill_list, ldestructive) + call util_spill(keeps%index2, discards%index2, lspill_list, ldestructive) + call util_spill(keeps%id1, discards%id1, lspill_list, ldestructive) + call util_spill(keeps%id2, discards%id2, lspill_list, ldestructive) + call util_spill(keeps%r1, discards%r1, lspill_list, ldestructive) + call util_spill(keeps%r2, discards%r2, lspill_list, ldestructive) + call util_spill(keeps%v1, discards%v1, lspill_list, ldestructive) + call util_spill(keeps%v2, discards%v2, lspill_list, ldestructive) + call util_spill(keeps%level, discards%level, lspill_list, ldestructive) nenc_old = keeps%nenc diff --git a/src/fraggle/fraggle_util.f90 b/src/fraggle/fraggle_util.f90 index 14a3b69df..3b0c42d6a 100644 --- a/src/fraggle/fraggle_util.f90 +++ b/src/fraggle/fraggle_util.f90 @@ -188,8 +188,8 @@ module subroutine fraggle_util_set_mass_dist(self, param) end if ! Sort the distribution in descending order by mass so that the largest fragment is always the first - call swiftest_util_sort(-mass, ind) - call swiftest_util_sort_rearrange(mass, ind, nfrag) + call util_sort(-mass, ind) + call util_sort_rearrange(mass, ind, nfrag) call move_alloc(mass, fragments%mass) fragments%Gmass(:) = G * fragments%mass(:) diff --git a/src/rmvs/rmvs_coarray.f90 b/src/rmvs/rmvs_coarray.f90 index 61c638c60..6f7467df8 100644 --- a/src/rmvs/rmvs_coarray.f90 +++ b/src/rmvs/rmvs_coarray.f90 @@ -31,6 +31,23 @@ module subroutine rmvs_coarray_coclone_cb(self) end subroutine rmvs_coarray_coclone_cb + module subroutine rmvs_coarray_coclone_interp(self) + !! author: David A. Minton + !! + !! Broadcasts the image 1 object to all other images in a coarray + implicit none + ! Arguments + class(rmvs_interp),intent(inout),codimension[*] :: self !! RMVS pl object + + call coclone(self%x) + call coclone(self%v) + call coclone(self%aobl) + call coclone(self%atide) + + return + end subroutine rmvs_coarray_coclone_interp + + module subroutine rmvs_coarray_coclone_pl(self) !! author: David A. Minton !! @@ -52,6 +69,26 @@ end subroutine rmvs_coarray_coclone_pl + module subroutine rmvs_coarray_coclone_system(self) + !! author: David A. Minton + !! + !! Broadcasts the image 1 object to all other images in a coarray + implicit none + ! Arguments + class(rmvs_nbody_system),intent(inout),codimension[*] :: self !! Swiftest body object + ! Internals + integer(I4B) :: i, img + + call coclone(self%lplanetocentric) + call coclone(self%rts) + call coclone(self%vbeg) + + call swiftest_coarray_coclone_system(self) + + return + end subroutine rmvs_coarray_coclone_system + + module subroutine rmvs_coarray_coclone_tp(self) !! author: David A. Minton !! @@ -84,10 +121,9 @@ module subroutine rmvs_coarray_component_clone_interp_arr1D(var,src_img) integer(I4B), intent(in),optional :: src_img ! Internals type(rmvs_interp), dimension(:), codimension[:], allocatable :: tmp - integer(I4B) :: img, si + integer(I4B) :: i,img, si integer(I4B), save :: n[*] logical, save :: isalloc[*] - if (present(src_img)) then si = src_img @@ -102,10 +138,14 @@ module subroutine rmvs_coarray_component_clone_interp_arr1D(var,src_img) if (.not. isalloc[si]) return allocate(tmp(n[si])[*]) + do i = 1, n[si] + call tmp(i)%coclone() + end do if (this_image() == si) then do img = 1, num_images() - tmp(:)[img] = var + tmp(:)[img] = var(:) end do + sync images(*) else sync images(si) @@ -128,8 +168,6 @@ module subroutine rmvs_coarray_cocollect_tp(self) call cocollect(self%lperi) call cocollect(self%plperP) call cocollect(self%plencP) - call cocollect(self%index) - call cocollect(self%ipleP) call swiftest_coarray_cocollect_tp(self) diff --git a/src/rmvs/rmvs_module.f90 b/src/rmvs/rmvs_module.f90 index 0f4c0f4b1..b219fd01b 100644 --- a/src/rmvs/rmvs_module.f90 +++ b/src/rmvs/rmvs_module.f90 @@ -35,6 +35,9 @@ module rmvs procedure :: dealloc => rmvs_util_dealloc_system !! Performs RMVS-specific deallocation procedure :: initialize => rmvs_util_setup_initialize_system !! Performs RMVS-specific initilization steps, including generating the close encounter planetocentric structures procedure :: step => rmvs_step_system !! Advance the RMVS nbody system forward in time by one step +#ifdef COARRAY + procedure :: coclone => rmvs_coarray_coclone_system !! Clones the image 1 body object to all other images in the coarray structure. +#endif end type rmvs_nbody_system type, private :: rmvs_interp @@ -45,6 +48,9 @@ module rmvs contains procedure :: dealloc => rmvs_util_dealloc_interp !! Deallocates all allocatable arrays final :: rmvs_final_interp !! Finalizes the RMVS interpolated nbody_system variables object - deallocates all allocatables +#ifdef COARRAY + procedure :: coclone => rmvs_coarray_coclone_interp +#endif end type rmvs_interp @@ -93,6 +99,7 @@ module rmvs final :: rmvs_final_tp !! Finalizes the RMVS test particle object - deallocates all allocatables #ifdef COARRAY procedure :: coclone => rmvs_coarray_coclone_tp !! Clones the image 1 body object to all other images in the coarray structure. + procedure :: cocollect => rmvs_coarray_cocollect_tp !! Clones the image 1 body object to all other images in the coarray structure. #endif end type rmvs_tp @@ -308,11 +315,21 @@ module subroutine rmvs_coarray_coclone_cb(self) class(rmvs_cb),intent(inout),codimension[*] :: self !! RMVS tp object end subroutine rmvs_coarray_coclone_cb + module subroutine rmvs_coarray_coclone_interp(self) + implicit none + class(rmvs_interp),intent(inout),codimension[*] :: self !! RMVS tp object + end subroutine rmvs_coarray_coclone_interp + module subroutine rmvs_coarray_coclone_pl(self) implicit none class(rmvs_pl),intent(inout),codimension[*] :: self !! RMVS pl object end subroutine rmvs_coarray_coclone_pl + module subroutine rmvs_coarray_coclone_system(self) + implicit none + class(rmvs_nbody_system),intent(inout),codimension[*] :: self !! RMVS nbody system object + end subroutine rmvs_coarray_coclone_system + module subroutine rmvs_coarray_coclone_tp(self) implicit none class(rmvs_tp),intent(inout),codimension[*] :: self !! RMVS tp object diff --git a/src/rmvs/rmvs_util.f90 b/src/rmvs/rmvs_util.f90 index 7100f7019..224ba2cd1 100644 --- a/src/rmvs/rmvs_util.f90 +++ b/src/rmvs/rmvs_util.f90 @@ -24,19 +24,17 @@ module subroutine rmvs_util_append_pl(self, source, lsource_mask) select type(source) class is (rmvs_pl) - associate(nold => self%nbody, nsrc => source%nbody) - call swiftest_util_append(self%nenc, source%nenc, nold, nsrc, lsource_mask) - call swiftest_util_append(self%tpenc1P, source%tpenc1P, nold, nsrc, lsource_mask) - call swiftest_util_append(self%plind, source%plind, nold, nsrc, lsource_mask) + call util_append(self%nenc, source%nenc, lsource_mask=lsource_mask) + call util_append(self%tpenc1P, source%tpenc1P, lsource_mask=lsource_mask) + call util_append(self%plind, source%plind, lsource_mask=lsource_mask) - ! The following are not implemented as RMVS doesn't make use of fill operations on pl type - ! So they are here as a placeholder in case someone wants to extend the RMVS class for some reason - !call swiftest_util_append(self%outer, source%outer, nold, nsrc, lsource_mask) - !call swiftest_util_append(self%inner, source%inner, nold, nsrc, lsource_mask) - !call swiftest_util_append(self%planetocentric, source%planetocentric, nold, nsrc, lsource_mask) + ! The following are not implemented as RMVS doesn't make use of fill operations on pl type + ! So they are here as a placeholder in case someone wants to extend the RMVS class for some reason + !call util_append(self%outer, source%outer, lsource_mask=lsource_mask) + !call util_append(self%inner, source%inner, lsource_mask=lsource_mask) + !call util_append(self%planetocentric, source%planetocentric, lsource_mask) - call whm_util_append_pl(self, source, lsource_mask) - end associate + call whm_util_append_pl(self, source, lsource_mask) class default write(*,*) "Invalid object passed to the append method. Source must be of class rmvs_pl or its descendents!" call base_util_exit(FAILURE) @@ -59,13 +57,11 @@ module subroutine rmvs_util_append_tp(self, source, lsource_mask) select type(source) class is (rmvs_tp) - associate(nold => self%nbody, nsrc => source%nbody) - call swiftest_util_append(self%lperi, source%lperi, nold, nsrc, lsource_mask) - call swiftest_util_append(self%plperP, source%plperP, nold, nsrc, lsource_mask) - call swiftest_util_append(self%plencP, source%plencP, nold, nsrc, lsource_mask) + call util_append(self%lperi, source%lperi, lsource_mask=lsource_mask) + call util_append(self%plperP, source%plperP, lsource_mask=lsource_mask) + call util_append(self%plencP, source%plencP, lsource_mask=lsource_mask) - call swiftest_util_append_tp(self, source, lsource_mask) ! Note: whm_tp does not have its own append method, so we skip back to the base class - end associate + call swiftest_util_append_tp(self, source, lsource_mask) ! Note: whm_tp does not have its own append method, so we skip back to the base class class default write(*,*) "Invalid object passed to the append method. Source must be of class rmvs_tp or its descendents!" call base_util_exit(FAILURE) @@ -177,15 +173,15 @@ module subroutine rmvs_util_fill_pl(self, inserts, lfill_list) associate(keeps => self) select type(inserts) class is (rmvs_pl) - call swiftest_util_fill(keeps%nenc, inserts%nenc, lfill_list) - call swiftest_util_fill(keeps%tpenc1P, inserts%tpenc1P, lfill_list) - call swiftest_util_fill(keeps%plind, inserts%plind, lfill_list) + call util_fill(keeps%nenc, inserts%nenc, lfill_list) + call util_fill(keeps%tpenc1P, inserts%tpenc1P, lfill_list) + call util_fill(keeps%plind, inserts%plind, lfill_list) ! The following are not implemented as RMVS doesn't make use of fill operations on pl type ! So they are here as a placeholder in case someone wants to extend the RMVS class for some reason - !call swiftest_util_fill(keeps%outer, inserts%outer, lfill_list) - !call swiftest_util_fill(keeps%inner, inserts%inner, lfill_list) - !call swiftest_util_fill(keeps%planetocentric, inserts%planetocentric, lfill_list) + !call util_fill(keeps%outer, inserts%outer, lfill_list) + !call util_fill(keeps%inner, inserts%inner, lfill_list) + !call util_fill(keeps%planetocentric, inserts%planetocentric, lfill_list) call whm_util_fill_pl(keeps, inserts, lfill_list) class default @@ -213,9 +209,9 @@ module subroutine rmvs_util_fill_tp(self, inserts, lfill_list) associate(keeps => self) select type(inserts) class is (rmvs_tp) - call swiftest_util_fill(keeps%lperi, inserts%lperi, lfill_list) - call swiftest_util_fill(keeps%plperP, inserts%plperP, lfill_list) - call swiftest_util_fill(keeps%plencP, inserts%plencP, lfill_list) + call util_fill(keeps%lperi, inserts%lperi, lfill_list) + call util_fill(keeps%plperP, inserts%plperP, lfill_list) + call util_fill(keeps%plencP, inserts%plencP, lfill_list) call swiftest_util_fill_tp(keeps, inserts, lfill_list) ! Note: whm_tp does not have its own fill method, so we skip back to the base class class default @@ -237,15 +233,15 @@ module subroutine rmvs_util_resize_pl(self, nnew) class(rmvs_pl), intent(inout) :: self !! RMVS massive body object integer(I4B), intent(in) :: nnew !! New size neded - call swiftest_util_resize(self%nenc, nnew) - call swiftest_util_resize(self%tpenc1P, nnew) - call swiftest_util_resize(self%plind, nnew) + call util_resize(self%nenc, nnew) + call util_resize(self%tpenc1P, nnew) + call util_resize(self%plind, nnew) ! The following are not implemented as RMVS doesn't make use of resize operations on pl type ! So they are here as a placeholder in case someone wants to extend the RMVS class for some reason - !call swiftest_util_resize(self%outer, nnew) - !call swiftest_util_resize(self%inner, nnew) - !call swiftest_util_resize(self%planetocentric, nnew) + !call util_resize(self%outer, nnew) + !call util_resize(self%inner, nnew) + !call util_resize(self%planetocentric, nnew) call whm_util_resize_pl(self, nnew) return @@ -261,10 +257,10 @@ module subroutine rmvs_util_resize_tp(self, nnew) class(rmvs_tp), intent(inout) :: self !! RMVS test particle object integer(I4B), intent(in) :: nnew !! New size neded - call swiftest_util_resize(self%lperi, nnew) - call swiftest_util_resize(self%plperP, nnew) - call swiftest_util_resize(self%plencP, nnew) - call swiftest_util_resize(self%rheliocentric, nnew) + call util_resize(self%lperi, nnew) + call util_resize(self%plperP, nnew) + call util_resize(self%plencP, nnew) + call util_resize(self%rheliocentric, nnew) call swiftest_util_resize_tp(self, nnew) @@ -453,11 +449,11 @@ module subroutine rmvs_util_sort_pl(self, sortby, ascending) associate(pl => self, npl => self%nbody) select case(sortby) case("nenc") - call swiftest_util_sort(direction * pl%nenc(1:npl), ind) + call util_sort(direction * pl%nenc(1:npl), ind) case("tpenc1P") - call swiftest_util_sort(direction * pl%tpenc1P(1:npl), ind) + call util_sort(direction * pl%tpenc1P(1:npl), ind) case("plind") - call swiftest_util_sort(direction * pl%plind(1:npl), ind) + call util_sort(direction * pl%plind(1:npl), ind) case("outer", "inner", "planetocentric", "lplanetocentric") write(*,*) 'Cannot sort by ' // trim(adjustl(sortby)) // '. Component not sortable!' case default ! Look for components in the parent class @@ -497,9 +493,9 @@ module subroutine rmvs_util_sort_tp(self, sortby, ascending) associate(tp => self, ntp => self%nbody) select case(sortby) case("plperP") - call swiftest_util_sort(direction * tp%plperP(1:ntp), ind) + call util_sort(direction * tp%plperP(1:ntp), ind) case("plencP") - call swiftest_util_sort(direction * tp%plencP(1:ntp), ind) + call util_sort(direction * tp%plencP(1:ntp), ind) case("lperi", "cb_heliocentric", "rheliocentric", "index", "ipleP", "lplanetocentric") write(*,*) 'Cannot sort by ' // trim(adjustl(sortby)) // '. Component not sortable!' case default ! Look for components in the parent class (*NOTE whm_tp does not need its own sort method, so we go straight to the swiftest_tp method) @@ -526,9 +522,9 @@ module subroutine rmvs_util_sort_rearrange_pl(self, ind) if (self%nbody == 0) return associate(pl => self, npl => self%nbody) - call swiftest_util_sort_rearrange(pl%nenc, ind, npl) - call swiftest_util_sort_rearrange(pl%tpenc1P, ind, npl) - call swiftest_util_sort_rearrange(pl%plind, ind, npl) + call util_sort_rearrange(pl%nenc, ind, npl) + call util_sort_rearrange(pl%tpenc1P, ind, npl) + call util_sort_rearrange(pl%plind, ind, npl) call swiftest_util_sort_rearrange_pl(pl,ind) end associate @@ -549,10 +545,10 @@ module subroutine rmvs_util_sort_rearrange_tp(self, ind) if (self%nbody == 0) return associate(tp => self, ntp => self%nbody) - call swiftest_util_sort_rearrange(tp%lperi, ind, ntp) - call swiftest_util_sort_rearrange(tp%plperP, ind, ntp) - call swiftest_util_sort_rearrange(tp%plencP, ind, ntp) - call swiftest_util_sort_rearrange(tp%rheliocentric, ind, ntp) + call util_sort_rearrange(tp%lperi, ind, ntp) + call util_sort_rearrange(tp%plperP, ind, ntp) + call util_sort_rearrange(tp%plencP, ind, ntp) + call util_sort_rearrange(tp%rheliocentric, ind, ntp) call swiftest_util_sort_rearrange_tp(tp,ind) end associate @@ -576,9 +572,9 @@ module subroutine rmvs_util_spill_pl(self, discards, lspill_list, ldestructive) associate(keeps => self) select type(discards) class is (rmvs_pl) - call swiftest_util_spill(keeps%nenc, discards%nenc, lspill_list, ldestructive) - call swiftest_util_spill(keeps%tpenc1P, discards%tpenc1P, lspill_list, ldestructive) - call swiftest_util_spill(keeps%plind, discards%plind, lspill_list, ldestructive) + call util_spill(keeps%nenc, discards%nenc, lspill_list, ldestructive) + call util_spill(keeps%tpenc1P, discards%tpenc1P, lspill_list, ldestructive) + call util_spill(keeps%plind, discards%plind, lspill_list, ldestructive) call whm_util_spill_pl(keeps, discards, lspill_list, ldestructive) class default @@ -607,9 +603,9 @@ module subroutine rmvs_util_spill_tp(self, discards, lspill_list, ldestructive) associate(keeps => self) select type(discards) class is (rmvs_tp) - call swiftest_util_spill(keeps%lperi, discards%lperi, lspill_list, ldestructive) - call swiftest_util_spill(keeps%plperP, discards%plperP, lspill_list, ldestructive) - call swiftest_util_spill(keeps%plencP, discards%plencP, lspill_list, ldestructive) + call util_spill(keeps%lperi, discards%lperi, lspill_list, ldestructive) + call util_spill(keeps%plperP, discards%plperP, lspill_list, ldestructive) + call util_spill(keeps%plencP, discards%plencP, lspill_list, ldestructive) call swiftest_util_spill_tp(keeps, discards, lspill_list, ldestructive) class default diff --git a/src/swiftest/swiftest_coarray.f90 b/src/swiftest/swiftest_coarray.f90 index 08590c112..2b04d0697 100644 --- a/src/swiftest/swiftest_coarray.f90 +++ b/src/swiftest/swiftest_coarray.f90 @@ -78,21 +78,19 @@ module subroutine swiftest_coarray_coclone_cb(self) call coclone(self%R0) call coclone(self%dR) - do i = 1, NDIM - call coclone(self%aobl(i)) - call coclone(self%atide(i)) - call coclone(self%aoblbeg(i)) - call coclone(self%aoblend(i)) - call coclone(self%atidebeg(i)) - call coclone(self%atideend(i)) - call coclone(self%rb(i)) - call coclone(self%vb(i)) - call coclone(self%agr(i)) - call coclone(self%Ip(i)) - call coclone(self%rot(i)) - call coclone(self%L0(i)) - call coclone(self%dL(i)) - end do + call coclonevec(self%aobl) + call coclonevec(self%atide) + call coclonevec(self%aoblbeg) + call coclonevec(self%aoblend) + call coclonevec(self%atidebeg) + call coclonevec(self%atideend) + call coclonevec(self%rb) + call coclonevec(self%vb) + call coclonevec(self%agr) + call coclonevec(self%Ip) + call coclonevec(self%rot) + call coclonevec(self%L0) + call coclonevec(self%dL) return end subroutine swiftest_coarray_coclone_cb @@ -120,8 +118,6 @@ module subroutine swiftest_coarray_coclone_pl(self) call coclone(self%k2) call coclone(self%Q ) call coclone(self%tlag) - call coclone(self%k_plpl) - call coclone(self%nplpl) call coclone(self%kin) call coclone(self%lmtiny) call coclone(self%nplm) @@ -143,8 +139,6 @@ module subroutine swiftest_coarray_coclone_tp(self) ! Arguments class(swiftest_tp),intent(inout),codimension[*] :: self !! Swiftest body object - call coclone(self%k_pltp) - call coclone(self%npltp) call coclone(self%nplenc) call swiftest_coarray_coclone_body(self) @@ -214,7 +208,8 @@ module subroutine swiftest_coarray_coclone_system(self) return end subroutine swiftest_coarray_coclone_system - + + module subroutine swiftest_coarray_component_clone_info(var,src_img) !! author: David A. Minton !! @@ -225,30 +220,31 @@ module subroutine swiftest_coarray_component_clone_info(var,src_img) type(swiftest_particle_info), intent(inout) :: var integer(I4B), intent(in),optional :: src_img ! Internals - type(swiftest_particle_info),save :: tmp[*] + type(swiftest_particle_info),allocatable :: tmp[:] integer(I4B) :: img, si + allocate(tmp[*]) if (present(src_img)) then - si = src_img + si = src_img else - si = 1 + si = 1 end if - + sync all if (this_image() == si) then - do img = 1, num_images() - tmp[img] = var - end do - sync images(*) + do img = 1, num_images() + tmp[img] = var + end do + sync images(*) else - sync images(si) - var = tmp[si] + sync images(si) + var = tmp[si] end if - + return end subroutine swiftest_coarray_component_clone_info - - + + module subroutine swiftest_coarray_component_clone_info_arr1D(var,src_img) !! author: David A. Minton !! @@ -261,33 +257,35 @@ module subroutine swiftest_coarray_component_clone_info_arr1D(var,src_img) ! Internals type(swiftest_particle_info), dimension(:), codimension[:], allocatable :: tmp integer(I4B) :: img, si - integer(I4B), save :: n[*] - logical, save :: isalloc[*] - + integer(I4B), allocatable :: n[:] + logical, allocatable :: isalloc[:] + + allocate(isalloc[*]) + allocate(n[*]) + if (present(src_img)) then - si = src_img + si = src_img else - si = 1 + si = 1 end if - - sync all + isalloc = allocated(var) if (isalloc) n = size(var) sync all if (.not. isalloc[si]) return - + allocate(tmp(n[si])[*]) if (this_image() == si) then - do img = 1, num_images() - tmp(:)[img] = var - end do - sync images(*) + do img = 1, num_images() + tmp(:)[img] = var + end do + sync images(*) else - sync images(si) - if (allocated(var)) deallocate(var) - allocate(var, source=tmp) + sync images(si) + if (allocated(var)) deallocate(var) + allocate(var, source=tmp) end if - + return end subroutine swiftest_coarray_component_clone_info_arr1D @@ -400,6 +398,8 @@ module subroutine swiftest_coarray_component_collect_info_arr1D(var,dest_img) return end subroutine swiftest_coarray_component_collect_info_arr1D + + module subroutine swiftest_coarray_cocollect_body(self) !! author: David A. Minton @@ -410,6 +410,9 @@ module subroutine swiftest_coarray_cocollect_body(self) class(swiftest_body),intent(inout), codimension[*] :: self !! Swiftest body object integer(I4B) :: i + if (this_image() == 1) write(*,*) "Before collect " + sync all + if (allocated(self%id)) write(*,*) "Image: ",this_image(), "id: ",self%id call cocollect(self%nbody) call cocollect(self%id) call cocollect(self%info) @@ -438,6 +441,10 @@ module subroutine swiftest_coarray_cocollect_body(self) call cocollect(self%omega) call cocollect(self%capm) + if (this_image() == 1) write(*,*) "after collect " + sync all + if (allocated(self%id)) write(*,*) "Image: ",this_image(), "id: ",self%id + return end subroutine swiftest_coarray_cocollect_body @@ -450,7 +457,6 @@ module subroutine swiftest_coarray_cocollect_tp(self) ! Arguments class(swiftest_tp),intent(inout),codimension[*] :: self !! Swiftest body object - call cocollect(self%k_pltp) call cocollect(self%npltp) call cocollect(self%nplenc) call swiftest_coarray_cocollect_body(self) @@ -488,9 +494,6 @@ module subroutine swiftest_coarray_collect_system(nbody_system, param) if (this_image() == 1) then write(param%display_unit,*) " Done collecting" - do i = 1, nbody_system%tp%nbody - write(*,*) i,"mu ",nbody_system%tp%mu(i) - end do end if return @@ -520,6 +523,10 @@ module subroutine swiftest_coarray_distribute_system(nbody_system, param) write(param%display_unit,*) " Distributing test particles across " // trim(adjustl(image_num_char)) // " images." end if + if (this_image() == 1) write(*,*) "Before distribute " + sync all + if (allocated(nbody_system%tp%id)) write(*,*) "Image: ",this_image(), "id: ",nbody_system%tp%id + ntp = nbody_system%tp%nbody sync all ntot = ntp[1] @@ -547,12 +554,15 @@ module subroutine swiftest_coarray_distribute_system(nbody_system, param) call nbody_system%tp%spill(tmp, lspill_list(:), ldestructive=.true.) deallocate(tmp, cotp) + + + if (this_image() == 1) write(*,*) "After distribute " + sync all + if (allocated(nbody_system%tp%id)) write(*,*) "Image: ",this_image(), "id: ",nbody_system%tp%id + if (this_image() == 1) then write(param%display_unit,*) " Done distributing" - do i = 1, nbody_system%tp%nbody - write(*,*) i,"mu ",nbody_system%tp%mu(i) - end do end if return @@ -584,7 +594,7 @@ module subroutine swiftest_coarray_initialize_system(nbody_system, param) end if return - end subroutine swiftest_coarray_initialize_system + end subroutine swiftest_coarray_initialize_system end submodule s_swiftest_coarray \ No newline at end of file diff --git a/src/swiftest/swiftest_driver.f90 b/src/swiftest/swiftest_driver.f90 index 9e139c986..4a34347f2 100644 --- a/src/swiftest/swiftest_driver.f90 +++ b/src/swiftest/swiftest_driver.f90 @@ -89,35 +89,47 @@ program swiftest_driver end if #endif -#ifdef COARRAY - if (this_image() == 1) then -#endif - call nbody_system%initialize(system_history, param) -#ifdef COARRAY - end if ! this_image() == 1 - call swiftest_coarray_initialize_system(nbody_system, param) +#ifdef COARRAY + ! The following line lets us read in the input files one image at a time + if (this_image() /= 1) sync images(this_image() - 1) +#endif + call nbody_system%initialize(system_history, param) +#ifdef COARRAY + if (this_image() == 1) then #endif ! If this is a new run, compute energy initial conditions (if energy tracking is turned on) and write the initial conditions to file. call nbody_system%display_run_information(param, integration_timer, phase="first") - if (param%lenergy) then - if (param%lrestart) then - call nbody_system%get_t0_values(system_history%nc, param) - else - call nbody_system%conservation_report(param, lterminal=.false.) ! This will save the initial values of energy and momentum - end if - call nbody_system%conservation_report(param, lterminal=.true.) + +#ifdef COARRAY + end if ! this_image() == 1 +#endif + if (param%lenergy) then + if (param%lrestart) then + call nbody_system%get_t0_values(system_history%nc, param) + else + call nbody_system%conservation_report(param, lterminal=.false.) ! This will save the initial values of energy and momentum end if - call system_history%take_snapshot(param,nbody_system) + call nbody_system%conservation_report(param, lterminal=.true.) + end if + call system_history%take_snapshot(param,nbody_system) + +#ifdef COARRAY + if (this_image() == 1) then +#endif call nbody_system%dump(param, system_history) - -#ifdef COARRAY +#ifdef COARRAY end if ! this_image() == 1 +#endif + +#ifdef COARRAY + if (this_image() < num_images()) sync images(this_image() + 1) ! Distribute test particles to the various images call nbody_system%coarray_distribute(param) #endif do iloop = istart, nloops !> Step the nbody_system forward in time + if (this_image() == 1) write(*,*) "Image: ", this_image(), "ntp: ",nbody_system%tp%nbody call integration_timer%start() call nbody_system%step(param, nbody_system%t, dt) call integration_timer%stop() diff --git a/src/swiftest/swiftest_io.f90 b/src/swiftest/swiftest_io.f90 index 71ba57b7f..68e60a04a 100644 --- a/src/swiftest/swiftest_io.f90 +++ b/src/swiftest/swiftest_io.f90 @@ -1556,7 +1556,7 @@ module subroutine swiftest_io_netcdf_write_frame_body(self, nc, param) associate(n => self%nbody, tslot => nc%tslot) if (n == 0) return - call swiftest_util_sort(self%id(1:n), ind) + call util_sort(self%id(1:n), ind) do i = 1, n j = ind(i) @@ -1750,7 +1750,7 @@ module subroutine swiftest_io_netcdf_write_info_body(self, nc, param) class is (swiftest_body) associate(n => self%nbody, tslot => nc%tslot) if (n == 0) return - call swiftest_util_sort(self%id(1:n), ind) + call util_sort(self%id(1:n), ind) call nc%get_idvals() do i = 1, n @@ -1871,7 +1871,7 @@ module subroutine swiftest_io_param_reader(self, unit, iotype, v_list, iostat, i integer(I4B), intent(out) :: iostat !! IO status code character(len=*), intent(inout) :: iomsg !! Message to pass if iostat /= 0 ! Internals - logical :: tstart_set = .false. !! Is the final time set in the input file? + logical :: tstart_set = .false. !! Is the final time set in the input file? logical :: tstop_set = .false. !! Is the final time set in the input file? logical :: dt_set = .false. !! Is the step size set in the input file? integer(I4B) :: ilength, ifirst, ilast, i !! Variables used to parse input file diff --git a/src/swiftest/swiftest_module.f90 b/src/swiftest/swiftest_module.f90 index 6388929ac..bbf1de296 100644 --- a/src/swiftest/swiftest_module.f90 +++ b/src/swiftest/swiftest_module.f90 @@ -1123,62 +1123,22 @@ module subroutine swiftest_user_kick_getacch_body(self, nbody_system, param, t, end subroutine swiftest_user_kick_getacch_body end interface - interface swiftest_util_append - module subroutine swiftest_util_append_arr_char_string(arr, source, nold, nsrc, lsource_mask) + interface util_append + module subroutine swiftest_util_append_arr_info(arr, source, nold, lsource_mask) implicit none - character(len=STRMAX), dimension(:), allocatable, intent(inout) :: arr !! Destination array - character(len=STRMAX), dimension(:), allocatable, intent(in) :: source !! Array to append - integer(I4B), intent(in) :: nold, nsrc !! Extend of the old array and the source array, respectively - logical, dimension(:), intent(in) :: lsource_mask !! Logical mask indicating which elements to append to - end subroutine swiftest_util_append_arr_char_string - - module subroutine swiftest_util_append_arr_DP(arr, source, nold, nsrc, lsource_mask) - implicit none - real(DP), dimension(:), allocatable, intent(inout) :: arr !! Destination array - real(DP), dimension(:), allocatable, intent(in) :: source !! Array to append - integer(I4B), intent(in) :: nold, nsrc !! Extend of the old array and the source array, respectively - logical, dimension(:), intent(in) :: lsource_mask !! Logical mask indicating which elements to append to - end subroutine swiftest_util_append_arr_DP - - module subroutine swiftest_util_append_arr_DPvec(arr, source, nold, nsrc, lsource_mask) - implicit none - real(DP), dimension(:,:), allocatable, intent(inout) :: arr !! Destination array - real(DP), dimension(:,:), allocatable, intent(in) :: source !! Array to append - integer(I4B), intent(in) :: nold, nsrc !! Extend of the old array and the source array, respectively - logical, dimension(:), intent(in) :: lsource_mask !! Logical mask indicating which elements to append to - end subroutine swiftest_util_append_arr_DPvec - - module subroutine swiftest_util_append_arr_I4B(arr, source, nold, nsrc, lsource_mask) - implicit none - integer(I4B), dimension(:), allocatable, intent(inout) :: arr !! Destination array - integer(I4B), dimension(:), allocatable, intent(in) :: source !! Array to append - integer(I4B), intent(in) :: nold, nsrc !! Extend of the old array and the source array, respectively - logical, dimension(:), intent(in) :: lsource_mask !! Logical mask indicating which elements to append to - end subroutine swiftest_util_append_arr_I4B - - module subroutine swiftest_util_append_arr_info(arr, source, nold, nsrc, lsource_mask) - implicit none - type(swiftest_particle_info), dimension(:), allocatable, intent(inout) :: arr !! Destination array - type(swiftest_particle_info), dimension(:), allocatable, intent(in) :: source !! Array to append - integer(I4B), intent(in) :: nold, nsrc !! Extend of the old array and the source array, respectively - logical, dimension(:), intent(in) :: lsource_mask !! Logical mask indicating which elements to append to + type(swiftest_particle_info), dimension(:), allocatable, intent(inout) :: arr !! Destination array + type(swiftest_particle_info), dimension(:), intent(in) :: source !! Array to append + integer(I4B), intent(in), optional :: nold !! Extent of original array. If passed, the source array will begin at arr(nold+1). Otherwise, the size of arr will be used. + logical, dimension(:), intent(in), optional :: lsource_mask !! Logical mask indicating which elements to append to end subroutine swiftest_util_append_arr_info - module subroutine swiftest_util_append_arr_kin(arr, source, nold, nsrc, lsource_mask) + module subroutine swiftest_util_append_arr_kin(arr, source, nold, lsource_mask) implicit none - type(swiftest_kinship), dimension(:), allocatable, intent(inout) :: arr !! Destination array - type(swiftest_kinship), dimension(:), allocatable, intent(in) :: source !! Array to append - integer(I4B), intent(in) :: nold, nsrc !! Extend of the old array and the source array, respectively - logical, dimension(:), intent(in) :: lsource_mask !! Logical mask indicating which elements to append to + type(swiftest_kinship), dimension(:), allocatable, intent(inout) :: arr !! Destination array + type(swiftest_kinship), dimension(:), intent(in) :: source !! Array to append + integer(I4B), intent(in), optional :: nold !! Extent of original array. If passed, the source array will begin at arr(nold+1). Otherwise, the size of arr will be used. + logical, dimension(:), intent(in), optional :: lsource_mask !! Logical mask indicating which elements to append to end subroutine swiftest_util_append_arr_kin - - module subroutine swiftest_util_append_arr_logical(arr, source, nold, nsrc, lsource_mask) - implicit none - logical, dimension(:), allocatable, intent(inout) :: arr !! Destination array - logical, dimension(:), allocatable, intent(in) :: source !! Array to append - integer(I4B), intent(in) :: nold, nsrc !! Extend of the old array and the source array, respectively - logical, dimension(:), intent(in) :: lsource_mask !! Logical mask indicating which elements to append to - end subroutine swiftest_util_append_arr_logical end interface interface @@ -1333,40 +1293,12 @@ module subroutine swiftest_util_fill_tp(self, inserts, lfill_list) end subroutine swiftest_util_fill_tp end interface - interface swiftest_util_fill - module subroutine swiftest_util_fill_arr_char_string(keeps, inserts, lfill_list) - implicit none - character(len=STRMAX), dimension(:), allocatable, intent(inout) :: keeps !! Array of values to keep - character(len=STRMAX), dimension(:), allocatable, intent(in) :: inserts !! Array of values to insert into keep - logical, dimension(:), intent(in) :: lfill_list !! Logical array of bodies to merge into the keeps - end subroutine swiftest_util_fill_arr_char_string - - module subroutine swiftest_util_fill_arr_DP(keeps, inserts, lfill_list) - implicit none - real(DP), dimension(:), allocatable, intent(inout) :: keeps !! Array of values to keep - real(DP), dimension(:), allocatable, intent(in) :: inserts !! Array of values to insert into keep - logical, dimension(:), intent(in) :: lfill_list !! Logical array of bodies to merge into the keeps - end subroutine swiftest_util_fill_arr_DP - - module subroutine swiftest_util_fill_arr_DPvec(keeps, inserts, lfill_list) - implicit none - real(DP), dimension(:,:), allocatable, intent(inout) :: keeps !! Array of values to keep - real(DP), dimension(:,:), allocatable, intent(in) :: inserts !! Array of values to insert into keep - logical, dimension(:), intent(in) :: lfill_list !! Logical array of bodies to merge into the keeps - end subroutine swiftest_util_fill_arr_DPvec - - module subroutine swiftest_util_fill_arr_I4B(keeps, inserts, lfill_list) - implicit none - integer(I4B), dimension(:), allocatable, intent(inout) :: keeps !! Array of values to keep - integer(I4B), dimension(:), allocatable, intent(in) :: inserts !! Array of values to insert into keep - logical, dimension(:), intent(in) :: lfill_list !! Logical array of bodies to merge into the keeps - end subroutine swiftest_util_fill_arr_I4B - + interface util_fill module subroutine swiftest_util_fill_arr_info(keeps, inserts, lfill_list) implicit none type(swiftest_particle_info), dimension(:), allocatable, intent(inout) :: keeps !! Array of values to keep type(swiftest_particle_info), dimension(:), allocatable, intent(in) :: inserts !! Array of values to insert into keep - logical, dimension(:), intent(in) :: lfill_list !! Logical array of bodies to merge into the keeps + logical, dimension(:), intent(in) :: lfill_list !! Logical array of bodies to merge into the keeps end subroutine swiftest_util_fill_arr_info module subroutine swiftest_util_fill_arr_kin(keeps, inserts, lfill_list) @@ -1375,13 +1307,6 @@ module subroutine swiftest_util_fill_arr_kin(keeps, inserts, lfill_list) type(swiftest_kinship), dimension(:), allocatable, intent(in) :: inserts !! Array of values to insert into keep logical, dimension(:), intent(in) :: lfill_list !! Logical array of bodies to merge into the keeps end subroutine swiftest_util_fill_arr_kin - - module subroutine swiftest_util_fill_arr_logical(keeps, inserts, lfill_list) - implicit none - logical, dimension(:), allocatable, intent(inout) :: keeps !! Array of values to keep - logical, dimension(:), allocatable, intent(in) :: inserts !! Array of values to insert into keep - logical, dimension(:), intent(in) :: lfill_list !! Logical array of bodies to merge into the keeps - end subroutine swiftest_util_fill_arr_logical end interface interface @@ -1526,31 +1451,7 @@ end subroutine swiftest_util_reset_kinship_pl end interface - interface swiftest_util_resize - module subroutine swiftest_util_resize_arr_char_string(arr, nnew) - implicit none - character(len=STRMAX), dimension(:), allocatable, intent(inout) :: arr !! Array to resize - integer(I4B), intent(in) :: nnew !! New size - end subroutine swiftest_util_resize_arr_char_string - - module subroutine swiftest_util_resize_arr_DP(arr, nnew) - implicit none - real(DP), dimension(:), allocatable, intent(inout) :: arr !! Array to resize - integer(I4B), intent(in) :: nnew !! New size - end subroutine swiftest_util_resize_arr_DP - - module subroutine swiftest_util_resize_arr_DPvec(arr, nnew) - implicit none - real(DP), dimension(:,:), allocatable, intent(inout) :: arr !! Array to resize - integer(I4B), intent(in) :: nnew !! New size - end subroutine swiftest_util_resize_arr_DPvec - - module subroutine swiftest_util_resize_arr_I4B(arr, nnew) - implicit none - integer(I4B), dimension(:), allocatable, intent(inout) :: arr !! Array to resize - integer(I4B), intent(in) :: nnew !! New size - end subroutine swiftest_util_resize_arr_I4B - + interface util_resize module subroutine swiftest_util_resize_arr_info(arr, nnew) implicit none type(swiftest_particle_info), dimension(:), allocatable, intent(inout) :: arr !! Array to resize @@ -1562,12 +1463,6 @@ module subroutine swiftest_util_resize_arr_kin(arr, nnew) type(swiftest_kinship), dimension(:), allocatable, intent(inout) :: arr !! Array to resize integer(I4B), intent(in) :: nnew !! New size end subroutine swiftest_util_resize_arr_kin - - module subroutine swiftest_util_resize_arr_logical(arr, nnew) - implicit none - logical, dimension(:), allocatable, intent(inout) :: arr !! Array to resize - integer(I4B), intent(in) :: nnew !! New size - end subroutine swiftest_util_resize_arr_logical end interface interface @@ -1686,89 +1581,8 @@ module subroutine swiftest_util_snapshot_system(self, param, nbody_system, t, ar end subroutine swiftest_util_snapshot_system end interface - interface swiftest_util_sort - pure module subroutine swiftest_util_sort_i4b(arr) - implicit none - integer(I4B), dimension(:), intent(inout) :: arr - end subroutine swiftest_util_sort_i4b - - pure module subroutine swiftest_util_sort_index_i4b(arr,ind) - implicit none - integer(I4B), dimension(:), intent(in) :: arr - integer(I4B), dimension(:), allocatable, intent(inout) :: ind - end subroutine swiftest_util_sort_index_i4b - - pure module subroutine swiftest_util_sort_index_I4B_I8Bind(arr,ind) - implicit none - integer(I4B), dimension(:), intent(in) :: arr - integer(I8B), dimension(:), allocatable, intent(inout) :: ind - end subroutine swiftest_util_sort_index_I4b_I8Bind - - pure module subroutine swiftest_util_sort_index_I8B_I8Bind(arr,ind) - implicit none - integer(I8B), dimension(:), intent(in) :: arr - integer(I8B), dimension(:), allocatable, intent(inout) :: ind - end subroutine swiftest_util_sort_index_I8B_I8Bind - - pure module subroutine swiftest_util_sort_sp(arr) - implicit none - real(SP), dimension(:), intent(inout) :: arr - end subroutine swiftest_util_sort_sp - - pure module subroutine swiftest_util_sort_index_sp(arr,ind) - implicit none - real(SP), dimension(:), intent(in) :: arr - integer(I4B), dimension(:), allocatable, intent(inout) :: ind - end subroutine swiftest_util_sort_index_sp - - pure module subroutine swiftest_util_sort_dp(arr) - implicit none - real(DP), dimension(:), intent(inout) :: arr - end subroutine swiftest_util_sort_dp - - pure module subroutine swiftest_util_sort_index_dp(arr,ind) - implicit none - real(DP), dimension(:), intent(in) :: arr - integer(I4B), dimension(:), allocatable, intent(inout) :: ind - end subroutine swiftest_util_sort_index_dp - end interface swiftest_util_sort - - interface swiftest_util_sort_rearrange - pure module subroutine swiftest_util_sort_rearrange_arr_char_string(arr, ind, n) - implicit none - character(len=STRMAX), dimension(:), allocatable, intent(inout) :: arr !! Destination array - integer(I4B), dimension(:), intent(in) :: ind !! Index to rearrange against - integer(I4B), intent(in) :: n !! Number of elements in arr and ind to rearrange - end subroutine swiftest_util_sort_rearrange_arr_char_string - - pure module subroutine swiftest_util_sort_rearrange_arr_DP(arr, ind, n) - implicit none - real(DP), dimension(:), allocatable, intent(inout) :: arr !! Destination array - integer(I4B), dimension(:), intent(in) :: ind !! Index to rearrange against - integer(I4B), intent(in) :: n !! Number of elements in arr and ind to rearrange - end subroutine swiftest_util_sort_rearrange_arr_DP - - pure module subroutine swiftest_util_sort_rearrange_arr_DPvec(arr, ind, n) - implicit none - real(DP), dimension(:,:), allocatable, intent(inout) :: arr !! Destination array - integer(I4B), dimension(:), intent(in) :: ind !! Index to rearrange against - integer(I4B), intent(in) :: n !! Number of elements in arr and ind to rearrange - end subroutine swiftest_util_sort_rearrange_arr_DPvec - - pure module subroutine swiftest_util_sort_rearrange_arr_I4B(arr, ind, n) - implicit none - integer(I4B), dimension(:), allocatable, intent(inout) :: arr !! Destination array - integer(I4B), dimension(:), intent(in) :: ind !! Index to rearrange against - integer(I4B), intent(in) :: n !! Number of elements in arr and ind to rearrange - end subroutine swiftest_util_sort_rearrange_arr_I4B - - pure module subroutine swiftest_util_sort_rearrange_arr_I4B_I8Bind(arr, ind, n) - implicit none - integer(I4B), dimension(:), allocatable, intent(inout) :: arr !! Destination array - integer(I8B), dimension(:), intent(in) :: ind !! Index to rearrange against - integer(I8B), intent(in) :: n !! Number of elements in arr and ind to rearrange - end subroutine swiftest_util_sort_rearrange_arr_I4B_I8Bind + interface util_sort_rearrange module subroutine swiftest_util_sort_rearrange_arr_info(arr, ind, n) implicit none type(swiftest_particle_info), dimension(:), allocatable, intent(inout) :: arr !! Destination array @@ -1783,20 +1597,7 @@ pure module subroutine swiftest_util_sort_rearrange_arr_kin(arr, ind, n) integer(I4B), intent(in) :: n !! Number of elements in arr and ind to rearrange end subroutine swiftest_util_sort_rearrange_arr_kin - pure module subroutine swiftest_util_sort_rearrange_arr_logical(arr, ind, n) - implicit none - logical, dimension(:), allocatable, intent(inout) :: arr !! Destination array - integer(I4B), dimension(:), intent(in) :: ind !! Index to rearrange against - integer(I4B), intent(in) :: n !! Number of elements in arr and ind to rearrange - end subroutine swiftest_util_sort_rearrange_arr_logical - - pure module subroutine swiftest_util_sort_rearrange_arr_logical_I8Bind(arr, ind, n) - implicit none - logical, dimension(:), allocatable, intent(inout) :: arr !! Destination array - integer(I8B), dimension(:), intent(in) :: ind !! Index to rearrange against - integer(I8B), intent(in) :: n !! Number of elements in arr and ind to rearrange - end subroutine swiftest_util_sort_rearrange_arr_logical_I8Bind - end interface swiftest_util_sort_rearrange + end interface util_sort_rearrange interface module subroutine swiftest_util_sort_rearrange_body(self, ind) @@ -1840,47 +1641,7 @@ end subroutine swiftest_util_sort_tp end interface - interface swiftest_util_spill - module subroutine swiftest_util_spill_arr_char_string(keeps, discards, lspill_list, ldestructive) - implicit none - character(len=STRMAX), dimension(:), allocatable, intent(inout) :: keeps !! Array of values to keep - character(len=STRMAX), dimension(:), allocatable, intent(inout) :: discards !! Array of discards - logical, dimension(:), intent(in) :: lspill_list !! Logical array of bodies to spill into the discardss - logical, intent(in) :: ldestructive !! Logical flag indicating whether or not this operation should alter the keeps array or not - end subroutine swiftest_util_spill_arr_char_string - - module subroutine swiftest_util_spill_arr_DP(keeps, discards, lspill_list, ldestructive) - implicit none - real(DP), dimension(:), allocatable, intent(inout) :: keeps !! Array of values to keep - real(DP), dimension(:), allocatable, intent(inout) :: discards !! Array of discards - logical, dimension(:), intent(in) :: lspill_list !! Logical array of bodies to spill into the discards - logical, intent(in) :: ldestructive !! Logical flag indicating whether or not this operation should alter the keeps array or not - end subroutine swiftest_util_spill_arr_DP - - module subroutine swiftest_util_spill_arr_DPvec(keeps, discards, lspill_list, ldestructive) - implicit none - real(DP), dimension(:,:), allocatable, intent(inout) :: keeps !! Array of values to keep - real(DP), dimension(:,:), allocatable, intent(inout) :: discards !! Array discards - logical, dimension(:), intent(in) :: lspill_list !! Logical array of bodies to spill into the discards - logical, intent(in) :: ldestructive !! Logical flag indicating whether or not this operation should alter the keeps array or not - end subroutine swiftest_util_spill_arr_DPvec - - module subroutine swiftest_util_spill_arr_I4B(keeps, discards, lspill_list, ldestructive) - implicit none - integer(I4B), dimension(:), allocatable, intent(inout) :: keeps !! Array of values to keep - integer(I4B), dimension(:), allocatable, intent(inout) :: discards !! Array of discards - logical, dimension(:), intent(in) :: lspill_list !! Logical array of bodies to spill into the discardss - logical, intent(in) :: ldestructive !! Logical flag indicating whether or not this operation should alter the keeps array or not - end subroutine swiftest_util_spill_arr_I4B - - module subroutine swiftest_util_spill_arr_I8B(keeps, discards, lspill_list, ldestructive) - implicit none - integer(I8B), dimension(:), allocatable, intent(inout) :: keeps !! Array of values to keep - integer(I8B), dimension(:), allocatable, intent(inout) :: discards !! Array of discards - logical, dimension(:), intent(in) :: lspill_list !! Logical array of bodies to spill into the discardss - logical, intent(in) :: ldestructive !! Logical flag indicating whether or not this operation should alter the keeps array or not - end subroutine swiftest_util_spill_arr_I8B - + interface util_spill module subroutine swiftest_util_spill_arr_info(keeps, discards, lspill_list, ldestructive) implicit none type(swiftest_particle_info), dimension(:), allocatable, intent(inout) :: keeps !! Array of values to keep @@ -1896,14 +1657,6 @@ module subroutine swiftest_util_spill_arr_kin(keeps, discards, lspill_list, ldes logical, dimension(:), intent(in) :: lspill_list !! Logical array of bodies to spill into the discardss logical, intent(in) :: ldestructive !! Logical flag indicating whether or not this operation should alter the keeps array or not end subroutine swiftest_util_spill_arr_kin - - module subroutine swiftest_util_spill_arr_logical(keeps, discards, lspill_list, ldestructive) - implicit none - logical, dimension(:), allocatable, intent(inout) :: keeps !! Array of values to keep - logical, dimension(:), allocatable, intent(inout) :: discards !! Array of discards - logical, dimension(:), intent(in) :: lspill_list !! Logical array of bodies to spill into the discardss - logical, intent(in) :: ldestructive !! Logical flag indicating whether or not this operation should alter the keeps array or not - end subroutine swiftest_util_spill_arr_logical end interface interface @@ -1933,22 +1686,6 @@ end subroutine swiftest_util_spill_tp end interface - interface swiftest_util_unique - module subroutine swiftest_util_unique_DP(input_array, output_array, index_map) - implicit none - real(DP), dimension(:), intent(in) :: input_array !! Unsorted input array - real(DP), dimension(:), allocatable, intent(out) :: output_array !! Sorted array of unique values - integer(I4B), dimension(:), allocatable, intent(out) :: index_map !! An array of the same size as input_array that such that any for any index i, output_array(index_map(i)) = input_array(i) - end subroutine swiftest_util_unique_DP - - module subroutine swiftest_util_unique_I4B(input_array, output_array, index_map) - implicit none - integer(I4B), dimension(:), intent(in) :: input_array !! Unsorted input array - integer(I4B), dimension(:), allocatable, intent(out) :: output_array !! Sorted array of unique values - integer(I4B), dimension(:), allocatable, intent(out) :: index_map !! An array of the same size as input_array that such that any for any index i, output_array(index_map(i)) = input_array(i) - end subroutine swiftest_util_unique_I4B - end interface swiftest_util_unique - interface module subroutine swiftest_util_valid_id_system(self, param) implicit none diff --git a/src/swiftest/swiftest_util.f90 b/src/swiftest/swiftest_util.f90 index 77747170b..fe46cc9e3 100644 --- a/src/swiftest/swiftest_util.f90 +++ b/src/swiftest/swiftest_util.f90 @@ -15,222 +15,97 @@ use fraggle contains - module subroutine swiftest_util_append_arr_char_string(arr, source, nold, nsrc, lsource_mask) + + module subroutine swiftest_util_append_arr_info(arr, source, nold, lsource_mask) !! author: David A. Minton !! - !! Append a single array of character string type onto another. If the destination array is not allocated, or is not big enough, this will allocate space for it. + !! Append a single array of particle information type onto another. If the destination array is not allocated, or is not big enough, this will allocate space for it. implicit none ! Arguments - character(len=STRMAX), dimension(:), allocatable, intent(inout) :: arr !! Destination array - character(len=STRMAX), dimension(:), allocatable, intent(in) :: source !! Array to append - integer(I4B), intent(in) :: nold, nsrc !! Extend of the old array and the source array, respectively - logical, dimension(:), intent(in) :: lsource_mask !! Logical mask indicating which elements to append to + type(swiftest_particle_info), dimension(:), allocatable, intent(inout) :: arr !! Destination array + type(swiftest_particle_info), dimension(:), intent(in) :: source !! Array to append + integer(I4B), intent(in), optional :: nold !! Extent of original array. If passed, the source array will begin at arr(nold+1). Otherwise, the size of arr will be used. + logical, dimension(:), intent(in), optional :: lsource_mask !! Logical mask indicating which elements to append to ! Internals - integer(I4B) :: nnew - - if (.not. allocated(source)) return - - nnew = count(lsource_mask(1:nsrc)) - if (nnew == 0) return + integer(I4B) :: nnew, nsrc, nend_orig, i + integer(I4B), dimension(:), allocatable :: idx - if (.not.allocated(arr)) then - allocate(arr(nold+nnew)) + if (present(lsource_mask)) then + nsrc = count(lsource_mask(:)) else - call swiftest_util_resize(arr, nold + nnew) + nsrc = size(source) end if - - arr(nold + 1:nold + nnew) = pack(source(1:nsrc), lsource_mask(1:nsrc)) - - return - end subroutine swiftest_util_append_arr_char_string - - - module subroutine swiftest_util_append_arr_DP(arr, source, nold, nsrc, lsource_mask) - !! author: David A. Minton - !! - !! Append a single array of double precision type onto another. If the destination array is not allocated, or is not big enough, this will allocate space for it. - implicit none - ! Arguments - real(DP), dimension(:), allocatable, intent(inout) :: arr !! Destination array - real(DP), dimension(:), allocatable, intent(in) :: source !! Array to append - integer(I4B), intent(in) :: nold, nsrc !! Extend of the old array and the source array, respectively - logical, dimension(:), intent(in) :: lsource_mask !! Logical mask indicating which elements to append to - ! Internals - integer(I4B) :: nnew - - if (.not. allocated(source)) return - - nnew = count(lsource_mask(1:nsrc)) - if (nnew == 0) return + if (nsrc == 0) return if (.not.allocated(arr)) then - allocate(arr(nold+nnew)) + nend_orig = 0 + allocate(arr(nsrc)) else - call swiftest_util_resize(arr, nold + nnew) + if (present(nold)) then + nend_orig = nold + else + nend_orig = size(arr) + end if + call util_resize(arr, nend_orig + nsrc) end if + nnew = nend_orig + nsrc - arr(nold + 1:nold + nnew) = pack(source(1:nsrc), lsource_mask(1:nsrc)) - - return - end subroutine swiftest_util_append_arr_DP - - - module subroutine swiftest_util_append_arr_DPvec(arr, source, nold, nsrc, lsource_mask) - !! author: David A. Minton - !! - !! Append a single array of double precision vector type of size (NDIM, n) onto another. If the destination array is not allocated, or is not big enough, this will allocate space for it. - implicit none - ! Arguments - real(DP), dimension(:,:), allocatable, intent(inout) :: arr !! Destination array - real(DP), dimension(:,:), allocatable, intent(in) :: source !! Array to append - integer(I4B), intent(in) :: nold, nsrc !! Extend of the old array and the source array, respectively - logical, dimension(:), intent(in) :: lsource_mask !! Logical mask indicating which elements to append to - ! Internals - integer(I4B) :: nnew - - if (.not. allocated(source)) return - - nnew = count(lsource_mask(1:nsrc)) - if (nnew == 0) return - - if (.not.allocated(arr)) then - allocate(arr(NDIM,nold+nnew)) + allocate(idx(nsrc)) + if (present(lsource_mask)) then + idx = pack([(i, i = 1, size(lsource_mask))], lsource_mask(:)) else - call swiftest_util_resize(arr, nold + nnew) - end if + idx = [(i, i = 1,nsrc)] + end if - arr(1, nold + 1:nold + nnew) = pack(source(1,1:nsrc), lsource_mask(1:nsrc)) - arr(2, nold + 1:nold + nnew) = pack(source(2,1:nsrc), lsource_mask(1:nsrc)) - arr(3, nold + 1:nold + nnew) = pack(source(3,1:nsrc), lsource_mask(1:nsrc)) + call swiftest_util_copy_particle_info_arr(source(:), arr(nold+1:nnew), idx) return - end subroutine swiftest_util_append_arr_DPvec + end subroutine swiftest_util_append_arr_info - module subroutine swiftest_util_append_arr_I4B(arr, source, nold, nsrc, lsource_mask) + module subroutine swiftest_util_append_arr_kin(arr, source, nold, lsource_mask) !! author: David A. Minton !! - !! Append a single array of integer(I4B) onto another. If the destination array is not allocated, or is not big enough, this will allocate space for it. + !! Append a single array of kinship type onto another. If the destination array is not allocated, or is not big enough, this will allocate space for it. implicit none ! Arguments - integer(I4B), dimension(:), allocatable, intent(inout) :: arr !! Destination array - integer(I4B), dimension(:), allocatable, intent(in) :: source !! Array to append - integer(I4B), intent(in) :: nold, nsrc !! Extend of the old array and the source array, respectively - logical, dimension(:), intent(in) :: lsource_mask !! Logical mask indicating which elements to append to + type(swiftest_kinship), dimension(:), allocatable, intent(inout) :: arr !! Destination array + type(swiftest_kinship), dimension(:), intent(in) :: source !! Array to append + integer(I4B), intent(in), optional :: nold !! Extent of original array. If passed, the source array will begin at arr(nold+1). Otherwise, the size of arr will be used. + logical, dimension(:), intent(in), optional :: lsource_mask !! Logical mask indicating which elements to append to ! Internals - integer(I4B) :: nnew + integer(I4B) :: nnew, nsrc, nend_orig - if (.not. allocated(source)) return - - nnew = count(lsource_mask(1:nsrc)) - if (nnew == 0) return - - if (.not.allocated(arr)) then - allocate(arr(nold+nnew)) + if (present(lsource_mask)) then + nsrc = count(lsource_mask(:)) else - call swiftest_util_resize(arr, nold + nnew) + nsrc = size(source) end if - - arr(nold + 1:nold + nnew) = pack(source(1:nsrc), lsource_mask(1:nsrc)) - - return - end subroutine swiftest_util_append_arr_I4B - - - module subroutine swiftest_util_append_arr_info(arr, source, nold, nsrc, lsource_mask) - !! author: David A. Minton - !! - !! Append a single array of particle information type onto another. If the destination array is not allocated, or is not big enough, this will allocate space for it. - implicit none - ! Arguments - type(swiftest_particle_info), dimension(:), allocatable, intent(inout) :: arr !! Destination array - type(swiftest_particle_info), dimension(:), allocatable, intent(in) :: source !! Array to append - integer(I4B), intent(in) :: nold, nsrc !! Extend of the old array and the source array, respectively - logical, dimension(:), intent(in) :: lsource_mask !! Logical mask indicating which elements to append to - ! Internals - integer(I4B) :: nnew, i - integer(I4B), dimension(:), allocatable :: idx - - if (.not. allocated(source)) return - - nnew = count(lsource_mask(1:nsrc)) - if (nnew == 0) return + if (nsrc == 0) return if (.not.allocated(arr)) then - allocate(arr(nold+nnew)) + nend_orig = 0 + allocate(arr(nsrc)) else - call swiftest_util_resize(arr, nold + nnew) + if (present(nold)) then + nend_orig = nold + else + nend_orig = size(arr) + end if + call util_resize(arr, nend_orig + nsrc) end if + nnew = nend_orig + nsrc - allocate(idx(nnew)) - - idx = pack([(i, i = 1, nsrc)], lsource_mask(1:nsrc)) - - call swiftest_util_copy_particle_info_arr(source(1:nsrc), arr(nold+1:nold+nnew), idx) - - return - end subroutine swiftest_util_append_arr_info - - - module subroutine swiftest_util_append_arr_kin(arr, source, nold, nsrc, lsource_mask) - !! author: David A. Minton - !! - !! Append a single array of kinship type onto another. If the destination array is not allocated, or is not big enough, this will allocate space for it. - implicit none - ! Arguments - type(swiftest_kinship), dimension(:), allocatable, intent(inout) :: arr !! Destination array - type(swiftest_kinship), dimension(:), allocatable, intent(in) :: source !! Array to append - integer(I4B), intent(in) :: nold, nsrc !! Extend of the old array and the source array, respectively - logical, dimension(:), intent(in) :: lsource_mask !! Logical mask indicating which elements to append to - ! Internals - integer(I4B) :: nnew - - if (.not. allocated(source)) return - - nnew = count(lsource_mask(1:nsrc)) - if (nnew == 0) return - - if (.not.allocated(arr)) then - allocate(arr(nold+nnew)) + if (present(lsource_mask)) then + arr(nold + 1:nnew) = pack(source(1:nsrc), lsource_mask(1:nsrc)) else - call swiftest_util_resize(arr, nold + nnew) + arr(nold + 1:nnew) = source(1:nsrc) end if - arr(nold + 1:nold + nnew) = pack(source(1:nsrc), lsource_mask(1:nsrc)) - return end subroutine swiftest_util_append_arr_kin - module subroutine swiftest_util_append_arr_logical(arr, source, nold, nsrc, lsource_mask) - !! author: David A. Minton - !! - !! Append a single array of logical type onto another. If the destination array is not allocated, or is not big enough, this will allocate space for it. - implicit none - ! Arguments - logical, dimension(:), allocatable, intent(inout) :: arr !! Destination array - logical, dimension(:), allocatable, intent(in) :: source !! Array to append - integer(I4B), intent(in) :: nold, nsrc !! Extend of the old array and the source array, respectively - logical, dimension(:), intent(in) :: lsource_mask !! Logical mask indicating which elements to append to - ! Internals - integer(I4B) :: nnew - - if (.not. allocated(source)) return - - nnew = count(lsource_mask(1:nsrc)) - if (nnew == 0) return - - if (.not.allocated(arr)) then - allocate(arr(nold+nnew)) - else - call swiftest_util_resize(arr, nold + nnew) - end if - - arr(nold + 1:nold + nnew) = pack(source(1:nsrc), lsource_mask(1:nsrc)) - - return - end subroutine swiftest_util_append_arr_logical - module subroutine swiftest_util_append_body(self, source, lsource_mask) !! author: David A. Minton @@ -243,40 +118,36 @@ module subroutine swiftest_util_append_body(self, source, lsource_mask) class(swiftest_body), intent(in) :: source !! Source object to append logical, dimension(:), intent(in) :: lsource_mask !! Logical mask indicating which elements to append to ! Internals - integer(I4B) :: nold, nsrc, nnew - - nold = self%nbody - nsrc = source%nbody - nnew = count(lsource_mask(1:nsrc)) - - call swiftest_util_append(self%id, source%id, nold, nsrc, lsource_mask) - call swiftest_util_append(self%info, source%info, nold, nsrc, lsource_mask) - call swiftest_util_append(self%lmask, source%lmask, nold, nsrc, lsource_mask) - call swiftest_util_append(self%status, source%status, nold, nsrc, lsource_mask) - call swiftest_util_append(self%ldiscard, source%ldiscard, nold, nsrc, lsource_mask) - call swiftest_util_append(self%lencounter, source%lencounter, nold, nsrc, lsource_mask) - call swiftest_util_append(self%lcollision, source%lcollision, nold, nsrc, lsource_mask) - call swiftest_util_append(self%mu, source%mu, nold, nsrc, lsource_mask) - call swiftest_util_append(self%rh, source%rh, nold, nsrc, lsource_mask) - call swiftest_util_append(self%vh, source%vh, nold, nsrc, lsource_mask) - call swiftest_util_append(self%rb, source%rb, nold, nsrc, lsource_mask) - call swiftest_util_append(self%vb, source%vb, nold, nsrc, lsource_mask) - call swiftest_util_append(self%ah, source%ah, nold, nsrc, lsource_mask) - call swiftest_util_append(self%aobl, source%aobl, nold, nsrc, lsource_mask) - call swiftest_util_append(self%atide, source%atide, nold, nsrc, lsource_mask) - call swiftest_util_append(self%agr, source%agr, nold, nsrc, lsource_mask) - call swiftest_util_append(self%ir3h, source%ir3h, nold, nsrc, lsource_mask) - call swiftest_util_append(self%isperi, source%isperi, nold, nsrc, lsource_mask) - call swiftest_util_append(self%peri, source%peri, nold, nsrc, lsource_mask) - call swiftest_util_append(self%atp, source%atp, nold, nsrc, lsource_mask) - call swiftest_util_append(self%a, source%a, nold, nsrc, lsource_mask) - call swiftest_util_append(self%e, source%e, nold, nsrc, lsource_mask) - call swiftest_util_append(self%inc, source%inc, nold, nsrc, lsource_mask) - call swiftest_util_append(self%capom, source%capom, nold, nsrc, lsource_mask) - call swiftest_util_append(self%omega, source%omega, nold, nsrc, lsource_mask) - call swiftest_util_append(self%capm, source%capm, nold, nsrc, lsource_mask) - - self%nbody = nold + nnew + + + call util_append(self%id, source%id, lsource_mask=lsource_mask) + call util_append(self%info, source%info, lsource_mask=lsource_mask) + call util_append(self%lmask, source%lmask, lsource_mask=lsource_mask) + call util_append(self%status, source%status, lsource_mask=lsource_mask) + call util_append(self%ldiscard, source%ldiscard, lsource_mask=lsource_mask) + call util_append(self%lencounter, source%lencounter, lsource_mask=lsource_mask) + call util_append(self%lcollision, source%lcollision, lsource_mask=lsource_mask) + call util_append(self%mu, source%mu, lsource_mask=lsource_mask) + call util_append(self%rh, source%rh, lsource_mask=lsource_mask) + call util_append(self%vh, source%vh, lsource_mask=lsource_mask) + call util_append(self%rb, source%rb, lsource_mask=lsource_mask) + call util_append(self%vb, source%vb, lsource_mask=lsource_mask) + call util_append(self%ah, source%ah, lsource_mask=lsource_mask) + call util_append(self%aobl, source%aobl, lsource_mask=lsource_mask) + call util_append(self%atide, source%atide, lsource_mask=lsource_mask) + call util_append(self%agr, source%agr, lsource_mask=lsource_mask) + call util_append(self%ir3h, source%ir3h, lsource_mask=lsource_mask) + call util_append(self%isperi, source%isperi, lsource_mask=lsource_mask) + call util_append(self%peri, source%peri, lsource_mask=lsource_mask) + call util_append(self%atp, source%atp, lsource_mask=lsource_mask) + call util_append(self%a, source%a, lsource_mask=lsource_mask) + call util_append(self%e, source%e, lsource_mask=lsource_mask) + call util_append(self%inc, source%inc, lsource_mask=lsource_mask) + call util_append(self%capom, source%capom, lsource_mask=lsource_mask) + call util_append(self%omega, source%omega, lsource_mask=lsource_mask) + call util_append(self%capm, source%capm, lsource_mask=lsource_mask) + + self%nbody = self%nbody + count(lsource_mask(:)) return end subroutine swiftest_util_append_body @@ -295,30 +166,28 @@ module subroutine swiftest_util_append_pl(self, source, lsource_mask) select type(source) class is (swiftest_pl) - associate(nold => self%nbody, nsrc => source%nbody) - call swiftest_util_append(self%mass, source%mass, nold, nsrc, lsource_mask) - call swiftest_util_append(self%Gmass, source%Gmass, nold, nsrc, lsource_mask) - call swiftest_util_append(self%rhill, source%rhill, nold, nsrc, lsource_mask) - call swiftest_util_append(self%renc, source%renc, nold, nsrc, lsource_mask) - call swiftest_util_append(self%radius, source%radius, nold, nsrc, lsource_mask) - call swiftest_util_append(self%density, source%density, nold, nsrc, lsource_mask) - call swiftest_util_append(self%rbeg, source%rbeg, nold, nsrc, lsource_mask) - call swiftest_util_append(self%rend, source%rend, nold, nsrc, lsource_mask) - call swiftest_util_append(self%vbeg, source%vbeg, nold, nsrc, lsource_mask) - call swiftest_util_append(self%Ip, source%Ip, nold, nsrc, lsource_mask) - call swiftest_util_append(self%rot, source%rot, nold, nsrc, lsource_mask) - call swiftest_util_append(self%k2, source%k2, nold, nsrc, lsource_mask) - call swiftest_util_append(self%Q, source%Q, nold, nsrc, lsource_mask) - call swiftest_util_append(self%tlag, source%tlag, nold, nsrc, lsource_mask) - call swiftest_util_append(self%kin, source%kin, nold, nsrc, lsource_mask) - call swiftest_util_append(self%lmtiny, source%lmtiny, nold, nsrc, lsource_mask) - call swiftest_util_append(self%nplenc, source%nplenc, nold, nsrc, lsource_mask) - call swiftest_util_append(self%ntpenc, source%ntpenc, nold, nsrc, lsource_mask) - - if (allocated(self%k_plpl)) deallocate(self%k_plpl) - - call swiftest_util_append_body(self, source, lsource_mask) - end associate + call util_append(self%mass, source%mass, lsource_mask=lsource_mask) + call util_append(self%Gmass, source%Gmass, lsource_mask=lsource_mask) + call util_append(self%rhill, source%rhill, lsource_mask=lsource_mask) + call util_append(self%renc, source%renc, lsource_mask=lsource_mask) + call util_append(self%radius, source%radius, lsource_mask=lsource_mask) + call util_append(self%density, source%density, lsource_mask=lsource_mask) + call util_append(self%rbeg, source%rbeg, lsource_mask=lsource_mask) + call util_append(self%rend, source%rend, lsource_mask=lsource_mask) + call util_append(self%vbeg, source%vbeg, lsource_mask=lsource_mask) + call util_append(self%Ip, source%Ip, lsource_mask=lsource_mask) + call util_append(self%rot, source%rot, lsource_mask=lsource_mask) + call util_append(self%k2, source%k2, lsource_mask=lsource_mask) + call util_append(self%Q, source%Q, lsource_mask=lsource_mask) + call util_append(self%tlag, source%tlag, lsource_mask=lsource_mask) + call util_append(self%kin, source%kin, lsource_mask=lsource_mask) + call util_append(self%lmtiny, source%lmtiny, lsource_mask=lsource_mask) + call util_append(self%nplenc, source%nplenc, lsource_mask=lsource_mask) + call util_append(self%ntpenc, source%ntpenc, lsource_mask=lsource_mask) + + if (allocated(self%k_plpl)) deallocate(self%k_plpl) + + call swiftest_util_append_body(self, source, lsource_mask) class default write(*,*) "Invalid object passed to the append method. Source must be of class swiftest_pl or its descendents" call base_util_exit(FAILURE) @@ -341,11 +210,9 @@ module subroutine swiftest_util_append_tp(self, source, lsource_mask) select type(source) class is (swiftest_tp) - associate(nold => self%nbody, nsrc => source%nbody) - call swiftest_util_append(self%nplenc, source%nplenc, nold, nsrc, lsource_mask) + call util_append(self%nplenc, source%nplenc, lsource_mask=lsource_mask) - call swiftest_util_append_body(self, source, lsource_mask) - end associate + call swiftest_util_append_body(self, source, lsource_mask) class default write(*,*) "Invalid object passed to the append method. Source must be of class swiftest_tp or its descendents" call base_util_exit(FAILURE) @@ -958,89 +825,6 @@ module subroutine swiftest_util_dealloc_tp(self) end subroutine swiftest_util_dealloc_tp - module subroutine swiftest_util_fill_arr_char_string(keeps, inserts, lfill_list) - !! author: David A. Minton - !! - !! Performs a fill operation on a single array of type character strings - !! This is the inverse of a spill operation - implicit none - ! Arguments - character(len=STRMAX), dimension(:), allocatable, intent(inout) :: keeps !! Array of values to keep - character(len=STRMAX), dimension(:), allocatable, intent(in) :: inserts !! Array of values to insert into keep - logical, dimension(:), intent(in) :: lfill_list !! Logical array of bodies to merge into the keeps - - if (.not.allocated(keeps) .or. .not.allocated(inserts)) return - - keeps(:) = unpack(keeps(:), .not.lfill_list(:), keeps(:)) - keeps(:) = unpack(inserts(:), lfill_list(:), keeps(:)) - - return - end subroutine swiftest_util_fill_arr_char_string - - - module subroutine swiftest_util_fill_arr_DP(keeps, inserts, lfill_list) - !! author: David A. Minton - !! - !! Performs a fill operation on a single array of type DP - !! This is the inverse of a spill operation - implicit none - ! Arguments - real(DP), dimension(:), allocatable, intent(inout) :: keeps !! Array of values to keep - real(DP), dimension(:), allocatable, intent(in) :: inserts !! Array of values to insert into keep - logical, dimension(:), intent(in) :: lfill_list !! Logical array of bodies to merge into the keeps - - if (.not.allocated(keeps) .or. .not.allocated(inserts)) return - - keeps(:) = unpack(keeps(:), .not.lfill_list(:), keeps(:)) - keeps(:) = unpack(inserts(:), lfill_list(:), keeps(:)) - - return - end subroutine swiftest_util_fill_arr_DP - - - module subroutine swiftest_util_fill_arr_DPvec(keeps, inserts, lfill_list) - !! author: David A. Minton - !! - !! Performs a fill operation on a single array of DP vectors with shape (NDIM, n) - !! This is the inverse of a spill operation - implicit none - ! Arguments - real(DP), dimension(:,:), allocatable, intent(inout) :: keeps !! Array of values to keep - real(DP), dimension(:,:), allocatable, intent(in) :: inserts !! Array of values to insert into keep - logical, dimension(:), intent(in) :: lfill_list !! Logical array of bodies to merge into the keeps - ! Internals - integer(I4B) :: i - - if (.not.allocated(keeps) .or. .not.allocated(inserts)) return - - do i = 1, NDIM - keeps(i,:) = unpack(keeps(i,:), .not.lfill_list(:), keeps(i,:)) - keeps(i,:) = unpack(inserts(i,:), lfill_list(:), keeps(i,:)) - end do - - return - end subroutine swiftest_util_fill_arr_DPvec - - - module subroutine swiftest_util_fill_arr_I4B(keeps, inserts, lfill_list) - !! author: David A. Minton - !! - !! Performs a fill operation on a single array of type I4B - !! This is the inverse of a spill operation - implicit none - ! Arguments - integer(I4B), dimension(:), allocatable, intent(inout) :: keeps !! Array of values to keep - integer(I4B), dimension(:), allocatable, intent(in) :: inserts !! Array of values to insert into keep - logical, dimension(:), intent(in) :: lfill_list !! Logical array of bodies to merge into the keeps - - if (.not.allocated(keeps) .or. .not.allocated(inserts)) return - - keeps(:) = unpack(keeps(:), .not.lfill_list(:), keeps(:)) - keeps(:) = unpack(inserts(:), lfill_list(:), keeps(:)) - - return - end subroutine swiftest_util_fill_arr_I4B - module subroutine swiftest_util_fill_arr_info(keeps, inserts, lfill_list) !! author: David A. Minton @@ -1070,26 +854,6 @@ module subroutine swiftest_util_fill_arr_info(keeps, inserts, lfill_list) end subroutine swiftest_util_fill_arr_info - module subroutine swiftest_util_fill_arr_logical(keeps, inserts, lfill_list) - !! author: David A. Minton - !! - !! Performs a fill operation on a single array of logicals - !! This is the inverse of a spill operation - implicit none - ! Arguments - logical, dimension(:), allocatable, intent(inout) :: keeps !! Array of values to keep - logical, dimension(:), allocatable, intent(in) :: inserts !! Array of values to insert into keep - logical, dimension(:), intent(in) :: lfill_list !! Logical array of bodies to merge into the keeps - - if (.not.allocated(keeps) .or. .not.allocated(inserts)) return - - keeps(:) = unpack(keeps(:), .not.lfill_list(:), keeps(:)) - keeps(:) = unpack(inserts(:), lfill_list(:), keeps(:)) - - return - end subroutine swiftest_util_fill_arr_logical - - module subroutine swiftest_util_fill_arr_kin(keeps, inserts, lfill_list) !! author: David A. Minton !! @@ -1124,32 +888,32 @@ module subroutine swiftest_util_fill_body(self, inserts, lfill_list) ! For each component, pack the discarded bodies into the discard object and do the inverse with the keeps !> Fill all the common components associate(keeps => self) - call swiftest_util_fill(keeps%id, inserts%id, lfill_list) - call swiftest_util_fill(keeps%info, inserts%info, lfill_list) - call swiftest_util_fill(keeps%lmask, inserts%lmask, lfill_list) - call swiftest_util_fill(keeps%status, inserts%status, lfill_list) - call swiftest_util_fill(keeps%ldiscard, inserts%ldiscard, lfill_list) - call swiftest_util_fill(keeps%lcollision, inserts%lcollision, lfill_list) - call swiftest_util_fill(keeps%lencounter, inserts%lencounter, lfill_list) - call swiftest_util_fill(keeps%mu, inserts%mu, lfill_list) - call swiftest_util_fill(keeps%rh, inserts%rh, lfill_list) - call swiftest_util_fill(keeps%vh, inserts%vh, lfill_list) - call swiftest_util_fill(keeps%rb, inserts%rb, lfill_list) - call swiftest_util_fill(keeps%vb, inserts%vb, lfill_list) - call swiftest_util_fill(keeps%ah, inserts%ah, lfill_list) - call swiftest_util_fill(keeps%aobl, inserts%aobl, lfill_list) - call swiftest_util_fill(keeps%agr, inserts%agr, lfill_list) - call swiftest_util_fill(keeps%atide, inserts%atide, lfill_list) - call swiftest_util_fill(keeps%ir3h, inserts%ir3h, lfill_list) - call swiftest_util_fill(keeps%isperi, inserts%isperi, lfill_list) - call swiftest_util_fill(keeps%peri, inserts%peri, lfill_list) - call swiftest_util_fill(keeps%atp, inserts%atp, lfill_list) - call swiftest_util_fill(keeps%a, inserts%a, lfill_list) - call swiftest_util_fill(keeps%e, inserts%e, lfill_list) - call swiftest_util_fill(keeps%inc, inserts%inc, lfill_list) - call swiftest_util_fill(keeps%capom, inserts%capom, lfill_list) - call swiftest_util_fill(keeps%omega, inserts%omega, lfill_list) - call swiftest_util_fill(keeps%capm, inserts%capm, lfill_list) + call util_fill(keeps%id, inserts%id, lfill_list) + call util_fill(keeps%info, inserts%info, lfill_list) + call util_fill(keeps%lmask, inserts%lmask, lfill_list) + call util_fill(keeps%status, inserts%status, lfill_list) + call util_fill(keeps%ldiscard, inserts%ldiscard, lfill_list) + call util_fill(keeps%lcollision, inserts%lcollision, lfill_list) + call util_fill(keeps%lencounter, inserts%lencounter, lfill_list) + call util_fill(keeps%mu, inserts%mu, lfill_list) + call util_fill(keeps%rh, inserts%rh, lfill_list) + call util_fill(keeps%vh, inserts%vh, lfill_list) + call util_fill(keeps%rb, inserts%rb, lfill_list) + call util_fill(keeps%vb, inserts%vb, lfill_list) + call util_fill(keeps%ah, inserts%ah, lfill_list) + call util_fill(keeps%aobl, inserts%aobl, lfill_list) + call util_fill(keeps%agr, inserts%agr, lfill_list) + call util_fill(keeps%atide, inserts%atide, lfill_list) + call util_fill(keeps%ir3h, inserts%ir3h, lfill_list) + call util_fill(keeps%isperi, inserts%isperi, lfill_list) + call util_fill(keeps%peri, inserts%peri, lfill_list) + call util_fill(keeps%atp, inserts%atp, lfill_list) + call util_fill(keeps%a, inserts%a, lfill_list) + call util_fill(keeps%e, inserts%e, lfill_list) + call util_fill(keeps%inc, inserts%inc, lfill_list) + call util_fill(keeps%capom, inserts%capom, lfill_list) + call util_fill(keeps%omega, inserts%omega, lfill_list) + call util_fill(keeps%capm, inserts%capm, lfill_list) ! This is the base class, so will be the last to be called in the cascade. keeps%nbody = size(keeps%id(:)) @@ -1175,23 +939,23 @@ module subroutine swiftest_util_fill_pl(self, inserts, lfill_list) select type (inserts) ! The standard requires us to select the type of both arguments in order to access all the components class is (swiftest_pl) !> Fill components specific to the massive body class - call swiftest_util_fill(keeps%mass, inserts%mass, lfill_list) - call swiftest_util_fill(keeps%Gmass, inserts%Gmass, lfill_list) - call swiftest_util_fill(keeps%rhill, inserts%rhill, lfill_list) - call swiftest_util_fill(keeps%renc, inserts%renc, lfill_list) - call swiftest_util_fill(keeps%radius, inserts%radius, lfill_list) - call swiftest_util_fill(keeps%density, inserts%density, lfill_list) - call swiftest_util_fill(keeps%rbeg, inserts%rbeg, lfill_list) - call swiftest_util_fill(keeps%rend, inserts%rend, lfill_list) - call swiftest_util_fill(keeps%vbeg, inserts%vbeg, lfill_list) - call swiftest_util_fill(keeps%Ip, inserts%Ip, lfill_list) - call swiftest_util_fill(keeps%rot, inserts%rot, lfill_list) - call swiftest_util_fill(keeps%k2, inserts%k2, lfill_list) - call swiftest_util_fill(keeps%Q, inserts%Q, lfill_list) - call swiftest_util_fill(keeps%tlag, inserts%tlag, lfill_list) - call swiftest_util_fill(keeps%kin, inserts%kin, lfill_list) - call swiftest_util_fill(keeps%nplenc, inserts%nplenc, lfill_list) - call swiftest_util_fill(keeps%ntpenc, inserts%ntpenc, lfill_list) + call util_fill(keeps%mass, inserts%mass, lfill_list) + call util_fill(keeps%Gmass, inserts%Gmass, lfill_list) + call util_fill(keeps%rhill, inserts%rhill, lfill_list) + call util_fill(keeps%renc, inserts%renc, lfill_list) + call util_fill(keeps%radius, inserts%radius, lfill_list) + call util_fill(keeps%density, inserts%density, lfill_list) + call util_fill(keeps%rbeg, inserts%rbeg, lfill_list) + call util_fill(keeps%rend, inserts%rend, lfill_list) + call util_fill(keeps%vbeg, inserts%vbeg, lfill_list) + call util_fill(keeps%Ip, inserts%Ip, lfill_list) + call util_fill(keeps%rot, inserts%rot, lfill_list) + call util_fill(keeps%k2, inserts%k2, lfill_list) + call util_fill(keeps%Q, inserts%Q, lfill_list) + call util_fill(keeps%tlag, inserts%tlag, lfill_list) + call util_fill(keeps%kin, inserts%kin, lfill_list) + call util_fill(keeps%nplenc, inserts%nplenc, lfill_list) + call util_fill(keeps%ntpenc, inserts%ntpenc, lfill_list) if (allocated(keeps%k_plpl)) deallocate(keeps%k_plpl) @@ -1220,7 +984,7 @@ module subroutine swiftest_util_fill_tp(self, inserts, lfill_list) select type(inserts) class is (swiftest_tp) !> Spill components specific to the test particle class - call swiftest_util_fill(keeps%nplenc, inserts%nplenc, lfill_list) + call util_fill(keeps%nplenc, inserts%nplenc, lfill_list) call swiftest_util_fill_body(keeps, inserts, lfill_list) class default @@ -1706,10 +1470,10 @@ module subroutine swiftest_util_index_map_storage(self) call swiftest_util_get_vals_storage(self, idvals, tvals) - call swiftest_util_unique(idvals,self%idvals,self%idmap) + call util_unique(idvals,self%idvals,self%idmap) self%nid = size(self%idvals) - call swiftest_util_unique(tvals,self%tvals,self%tmap) + call util_unique(tvals,self%tvals,self%tmap) self%nt = size(self%tvals) return @@ -2044,16 +1808,16 @@ module subroutine swiftest_util_reset_kinship_pl(self, idx) end subroutine swiftest_util_reset_kinship_pl - module subroutine swiftest_util_resize_arr_char_string(arr, nnew) + module subroutine swiftest_util_resize_arr_info(arr, nnew) !! author: David A. Minton !! - !! Resizes an array component of type character string. nnew = 0 will deallocate. + !! Resizes an array component of type character string. Array will only be resized if has previously been allocated. Passing nnew = 0 will deallocate. implicit none ! Arguments - character(len=STRMAX), dimension(:), allocatable, intent(inout) :: arr !! Array to resize - integer(I4B), intent(in) :: nnew !! New size + type(swiftest_particle_info), dimension(:), allocatable, intent(inout) :: arr !! Array to resize + integer(I4B), intent(in) :: nnew !! New size ! Internals - character(len=STRMAX), dimension(:), allocatable :: tmp !! Temporary storage array in case the input array is already allocated + type(swiftest_particle_info), dimension(:), allocatable :: tmp !! Temporary storage array in case the input array is already allocated integer(I4B) :: nold !! Old size if (nnew < 0) return @@ -2072,34 +1836,29 @@ module subroutine swiftest_util_resize_arr_char_string(arr, nnew) if (nnew == nold) return allocate(tmp(nnew)) - if (nold > 0) then - if (nnew > nold) then - tmp(1:nold) = arr(1:nold) - tmp(nold+1:nnew) = "" - else - tmp(1:nnew) = arr(1:nnew) - end if + if (nnew > nold) then + call swiftest_util_copy_particle_info_arr(arr(1:nold), tmp(1:nold)) else - tmp(1:nnew) = "" + call swiftest_util_copy_particle_info_arr(arr(1:nnew), tmp(1:nnew)) end if + call move_alloc(tmp, arr) return - end subroutine swiftest_util_resize_arr_char_string - + end subroutine swiftest_util_resize_arr_info - module subroutine swiftest_util_resize_arr_DP(arr, nnew) + + module subroutine swiftest_util_resize_arr_kin(arr, nnew) !! author: David A. Minton !! - !! Resizes an array component of double precision type. Passing nnew = 0 will deallocate. + !! Resizes an array component of type character string. Array will only be resized if has previously been allocated. Passing nnew = 0 will deallocate. implicit none ! Arguments - real(DP), dimension(:), allocatable, intent(inout) :: arr !! Array to resize - integer(I4B), intent(in) :: nnew !! New size + type(swiftest_kinship), dimension(:), allocatable, intent(inout) :: arr !! Array to resize + integer(I4B), intent(in) :: nnew !! New size ! Internals - real(DP), dimension(:), allocatable :: tmp !! Temporary storage array in case the input array is already allocated + type(swiftest_kinship), dimension(:), allocatable :: tmp !! Temporary storage array in case the input array is already allocated integer(I4B) :: nold !! Old size - real(DP), parameter :: init_val = 0.0_DP if (nnew < 0) return @@ -2114,324 +1873,98 @@ module subroutine swiftest_util_resize_arr_DP(arr, nnew) nold = 0 end if - if (nnew == nold) return - allocate(tmp(nnew)) - if (nold > 0) then - if (nnew > nold) then - tmp(1:nold) = arr(1:nold) - tmp(nold+1:nnew) = init_val - else - tmp(1:nnew) = arr(1:nnew) - end if + if (nnew > nold) then + tmp(1:nold) = arr(1:nold) else - tmp(1:nnew) = init_val + tmp(1:nnew) = arr(1:nnew) end if call move_alloc(tmp, arr) return - end subroutine swiftest_util_resize_arr_DP + end subroutine swiftest_util_resize_arr_kin - module subroutine swiftest_util_resize_arr_DPvec(arr, nnew) + module subroutine swiftest_util_resize_body(self, nnew) !! author: David A. Minton !! - !! Resizes an array component of double precision vectors of size (NDIM, n). Passing nnew = 0 will deallocate. + !! Checks the current size of a Swiftest body against the requested size and resizes it if it is too small. implicit none ! Arguments - real(DP), dimension(:,:), allocatable, intent(inout) :: arr !! Array to resize - integer(I4B), intent(in) :: nnew !! New size - ! Internals - real(DP), dimension(:,:), allocatable :: tmp !! Temporary storage array in case the input array is already allocated - integer(I4B) :: nold !! Old size - real(DP), dimension(NDIM), parameter :: init_val = 0.0_DP - integer(I4B) :: i - - if (nnew < 0) return - - if (nnew == 0) then - if (allocated(arr)) deallocate(arr) - return - end if - - if (allocated(arr)) then - nold = size(arr, dim=2) - else - nold = 0 - end if - - if (nnew == nold) return - - allocate(tmp(NDIM, nnew)) - if (nold > 0) then - if (nnew > nold) then - tmp(:,1:nold) = arr(:,1:nold) - do i = nold+1, nnew - tmp(:,i) = init_val(:) - end do - else - tmp(:,1:nnew) = arr(:,1:nnew) - end if - else - do i = 1, nnew - tmp(:, i) = init_val(:) - end do - end if - call move_alloc(tmp, arr) + class(swiftest_body), intent(inout) :: self !! Swiftest body object + integer(I4B), intent(in) :: nnew !! New size neded - return + call util_resize(self%info, nnew) + call util_resize(self%id, nnew) + call util_resize(self%status, nnew) + call util_resize(self%lcollision, nnew) + call util_resize(self%lencounter, nnew) + call util_resize(self%ldiscard, nnew) + call util_resize(self%lmask, nnew) + call util_resize(self%mu, nnew) + call util_resize(self%rh, nnew) + call util_resize(self%vh, nnew) + call util_resize(self%rb, nnew) + call util_resize(self%vb, nnew) + call util_resize(self%ah, nnew) + call util_resize(self%aobl, nnew) + call util_resize(self%atide, nnew) + call util_resize(self%agr, nnew) + call util_resize(self%ir3h, nnew) + call util_resize(self%a, nnew) + call util_resize(self%e, nnew) + call util_resize(self%inc, nnew) + call util_resize(self%capom, nnew) + call util_resize(self%omega, nnew) + call util_resize(self%capm, nnew) + self%nbody = count(self%status(1:nnew) /= INACTIVE) return - end subroutine swiftest_util_resize_arr_DPvec + end subroutine swiftest_util_resize_body - module subroutine swiftest_util_resize_arr_I4B(arr, nnew) + module subroutine swiftest_util_resize_pl(self, nnew) !! author: David A. Minton !! - !! Resizes an array component of integer type. Passing nnew = 0 will deallocate. + !! Checks the current size of a Swiftest massive body against the requested size and resizes it if it is too small. implicit none ! Arguments - integer(I4B), dimension(:), allocatable, intent(inout) :: arr !! Array to resize - integer(I4B), intent(in) :: nnew !! New size - ! Internals - integer(I4B), dimension(:), allocatable :: tmp !! Temporary storage array in case the input array is already allocated - integer(I4B) :: nold !! Old size - integer(I4B), parameter :: init_val = -1 + class(swiftest_pl), intent(inout) :: self !! Swiftest massive body object + integer(I4B), intent(in) :: nnew !! New size neded - if (nnew < 0) return + call swiftest_util_resize_body(self, nnew) - if (nnew == 0) then - if (allocated(arr)) deallocate(arr) - return - end if - - if (allocated(arr)) then - nold = size(arr) - else - nold = 0 - end if + call util_resize(self%mass, nnew) + call util_resize(self%Gmass, nnew) + call util_resize(self%rhill, nnew) + call util_resize(self%renc, nnew) + call util_resize(self%radius, nnew) + call util_resize(self%rbeg, nnew) + call util_resize(self%rend, nnew) + call util_resize(self%vbeg, nnew) + call util_resize(self%density, nnew) + call util_resize(self%Ip, nnew) + call util_resize(self%rot, nnew) + call util_resize(self%k2, nnew) + call util_resize(self%Q, nnew) + call util_resize(self%tlag, nnew) + call util_resize(self%kin, nnew) + call util_resize(self%lmtiny, nnew) + call util_resize(self%nplenc, nnew) + call util_resize(self%ntpenc, nnew) - if (nnew == nold) return - - allocate(tmp(nnew)) - if (nold > 0) then - if (nnew > nold) then - tmp(1:nold) = arr(1:nold) - tmp(nold+1:nnew) = init_val - else - tmp(1:nnew) = arr(1:nnew) - end if - else - tmp(1:nnew) = init_val - end if - call move_alloc(tmp, arr) + + + if (allocated(self%k_plpl)) deallocate(self%k_plpl) return - end subroutine swiftest_util_resize_arr_I4B + end subroutine swiftest_util_resize_pl - module subroutine swiftest_util_resize_arr_info(arr, nnew) + module subroutine swiftest_util_resize_tp(self, nnew) !! author: David A. Minton !! - !! Resizes an array component of type character string. Array will only be resized if has previously been allocated. Passing nnew = 0 will deallocate. - implicit none - ! Arguments - type(swiftest_particle_info), dimension(:), allocatable, intent(inout) :: arr !! Array to resize - integer(I4B), intent(in) :: nnew !! New size - ! Internals - type(swiftest_particle_info), dimension(:), allocatable :: tmp !! Temporary storage array in case the input array is already allocated - integer(I4B) :: nold !! Old size - - if (nnew < 0) return - - if (nnew == 0) then - if (allocated(arr)) deallocate(arr) - return - end if - - if (allocated(arr)) then - nold = size(arr) - else - nold = 0 - end if - - if (nnew == nold) return - - allocate(tmp(nnew)) - if (nnew > nold) then - call swiftest_util_copy_particle_info_arr(arr(1:nold), tmp(1:nold)) - else - call swiftest_util_copy_particle_info_arr(arr(1:nnew), tmp(1:nnew)) - end if - - call move_alloc(tmp, arr) - - return - end subroutine swiftest_util_resize_arr_info - - - module subroutine swiftest_util_resize_arr_kin(arr, nnew) - !! author: David A. Minton - !! - !! Resizes an array component of type character string. Array will only be resized if has previously been allocated. Passing nnew = 0 will deallocate. - implicit none - ! Arguments - type(swiftest_kinship), dimension(:), allocatable, intent(inout) :: arr !! Array to resize - integer(I4B), intent(in) :: nnew !! New size - ! Internals - type(swiftest_kinship), dimension(:), allocatable :: tmp !! Temporary storage array in case the input array is already allocated - integer(I4B) :: nold !! Old size - - if (nnew < 0) return - - if (nnew == 0) then - if (allocated(arr)) deallocate(arr) - return - end if - - if (allocated(arr)) then - nold = size(arr) - else - nold = 0 - end if - - allocate(tmp(nnew)) - if (nnew > nold) then - tmp(1:nold) = arr(1:nold) - else - tmp(1:nnew) = arr(1:nnew) - end if - call move_alloc(tmp, arr) - - return - end subroutine swiftest_util_resize_arr_kin - - - module subroutine swiftest_util_resize_arr_logical(arr, nnew) - !! author: David A. Minton - !! - !! Resizes an array component of logical type. Passing nnew = 0 will deallocate. - implicit none - ! Arguments - logical, dimension(:), allocatable, intent(inout) :: arr !! Array to resize - integer(I4B), intent(in) :: nnew !! New size - ! Internals - logical, dimension(:), allocatable :: tmp !! Temporary storage array in case the input array is already allocated - integer(I4B) :: nold !! Old size - logical, parameter :: init_val = .false. - - if (nnew < 0) return - - if (nnew == 0) then - if (allocated(arr)) deallocate(arr) - return - end if - - if (allocated(arr)) then - nold = size(arr) - else - nold = 0 - end if - - if (nnew == nold) return - - allocate(tmp(nnew)) - if (nold > 0) then - if (nnew > nold) then - tmp(1:nold) = arr(1:nold) - tmp(nold+1:nnew) = init_val - else - tmp(1:nnew) = arr(1:nnew) - end if - else - tmp(1:nnew) = init_val - end if - call move_alloc(tmp, arr) - - return - end subroutine swiftest_util_resize_arr_logical - - - module subroutine swiftest_util_resize_body(self, nnew) - !! author: David A. Minton - !! - !! Checks the current size of a Swiftest body against the requested size and resizes it if it is too small. - implicit none - ! Arguments - class(swiftest_body), intent(inout) :: self !! Swiftest body object - integer(I4B), intent(in) :: nnew !! New size neded - - call swiftest_util_resize(self%info, nnew) - call swiftest_util_resize(self%id, nnew) - call swiftest_util_resize(self%status, nnew) - call swiftest_util_resize(self%lcollision, nnew) - call swiftest_util_resize(self%lencounter, nnew) - call swiftest_util_resize(self%ldiscard, nnew) - call swiftest_util_resize(self%lmask, nnew) - call swiftest_util_resize(self%mu, nnew) - call swiftest_util_resize(self%rh, nnew) - call swiftest_util_resize(self%vh, nnew) - call swiftest_util_resize(self%rb, nnew) - call swiftest_util_resize(self%vb, nnew) - call swiftest_util_resize(self%ah, nnew) - call swiftest_util_resize(self%aobl, nnew) - call swiftest_util_resize(self%atide, nnew) - call swiftest_util_resize(self%agr, nnew) - call swiftest_util_resize(self%ir3h, nnew) - call swiftest_util_resize(self%a, nnew) - call swiftest_util_resize(self%e, nnew) - call swiftest_util_resize(self%inc, nnew) - call swiftest_util_resize(self%capom, nnew) - call swiftest_util_resize(self%omega, nnew) - call swiftest_util_resize(self%capm, nnew) - self%nbody = count(self%status(1:nnew) /= INACTIVE) - - return - end subroutine swiftest_util_resize_body - - - module subroutine swiftest_util_resize_pl(self, nnew) - !! author: David A. Minton - !! - !! Checks the current size of a Swiftest massive body against the requested size and resizes it if it is too small. - implicit none - ! Arguments - class(swiftest_pl), intent(inout) :: self !! Swiftest massive body object - integer(I4B), intent(in) :: nnew !! New size neded - - call swiftest_util_resize_body(self, nnew) - - call swiftest_util_resize(self%mass, nnew) - call swiftest_util_resize(self%Gmass, nnew) - call swiftest_util_resize(self%rhill, nnew) - call swiftest_util_resize(self%renc, nnew) - call swiftest_util_resize(self%radius, nnew) - call swiftest_util_resize(self%rbeg, nnew) - call swiftest_util_resize(self%rend, nnew) - call swiftest_util_resize(self%vbeg, nnew) - call swiftest_util_resize(self%density, nnew) - call swiftest_util_resize(self%Ip, nnew) - call swiftest_util_resize(self%rot, nnew) - call swiftest_util_resize(self%k2, nnew) - call swiftest_util_resize(self%Q, nnew) - call swiftest_util_resize(self%tlag, nnew) - call swiftest_util_resize(self%kin, nnew) - call swiftest_util_resize(self%lmtiny, nnew) - call swiftest_util_resize(self%nplenc, nnew) - call swiftest_util_resize(self%ntpenc, nnew) - - - - if (allocated(self%k_plpl)) deallocate(self%k_plpl) - - return - end subroutine swiftest_util_resize_pl - - - module subroutine swiftest_util_resize_tp(self, nnew) - !! author: David A. Minton - !! - !! Checks the current size of a Swiftest test particle against the requested size and resizes it if it is too small. + !! Checks the current size of a Swiftest test particle against the requested size and resizes it if it is too small. implicit none ! Arguments class(swiftest_tp), intent(inout) :: self !! Swiftest test particle object @@ -2439,10 +1972,10 @@ module subroutine swiftest_util_resize_tp(self, nnew) call swiftest_util_resize_body(self, nnew) - call swiftest_util_resize(self%nplenc, nnew) - call swiftest_util_resize(self%isperi, nnew) - call swiftest_util_resize(self%peri, nnew) - call swiftest_util_resize(self%atp, nnew) + call util_resize(self%nplenc, nnew) + call util_resize(self%isperi, nnew) + call util_resize(self%peri, nnew) + call util_resize(self%atp, nnew) return end subroutine swiftest_util_resize_tp @@ -3154,25 +2687,25 @@ module subroutine swiftest_util_sort_body(self, sortby, ascending) associate(body => self, n => self%nbody) select case(sortby) case("id") - call swiftest_util_sort(direction * body%id(1:n), ind) + call util_sort(direction * body%id(1:n), ind) case("status") - call swiftest_util_sort(direction * body%status(1:n), ind) + call util_sort(direction * body%status(1:n), ind) case("ir3h") - call swiftest_util_sort(direction * body%ir3h(1:n), ind) + call util_sort(direction * body%ir3h(1:n), ind) case("a") - call swiftest_util_sort(direction * body%a(1:n), ind) + call util_sort(direction * body%a(1:n), ind) case("e") - call swiftest_util_sort(direction * body%e(1:n), ind) + call util_sort(direction * body%e(1:n), ind) case("inc") - call swiftest_util_sort(direction * body%inc(1:n), ind) + call util_sort(direction * body%inc(1:n), ind) case("capom") - call swiftest_util_sort(direction * body%capom(1:n), ind) + call util_sort(direction * body%capom(1:n), ind) case("mu") - call swiftest_util_sort(direction * body%mu(1:n), ind) + call util_sort(direction * body%mu(1:n), ind) case("peri") - call swiftest_util_sort(direction * body%peri(1:n), ind) + call util_sort(direction * body%peri(1:n), ind) case("atp") - call swiftest_util_sort(direction * body%atp(1:n), ind) + call util_sort(direction * body%atp(1:n), ind) case("info", "lfirst", "nbody", "ldiscard", "lcollision", "lencounter", "rh", "vh", "rb", "vb", "ah", "aobl", "atide", "agr","isperi") write(*,*) 'Cannot sort by ' // trim(adjustl(sortby)) // '. Component not sortable!' case default @@ -3188,856 +2721,167 @@ module subroutine swiftest_util_sort_body(self, sortby, ascending) end subroutine swiftest_util_sort_body - pure module subroutine swiftest_util_sort_dp(arr) + module subroutine swiftest_util_sort_pl(self, sortby, ascending) !! author: David A. Minton !! - !! Sort input DP precision array in place into ascending numerical order using quicksort. - !! + !! Sort a Swiftest massive body object in-place. + !! sortby is a string indicating which array component to sort. implicit none ! Arguments - real(DP), dimension(:), intent(inout) :: arr + class(swiftest_pl), intent(inout) :: self !! Swiftest massive body object + character(*), intent(in) :: sortby !! Sorting attribute + logical, intent(in) :: ascending !! Logical flag indicating whether or not the sorting should be in ascending or descending order + ! Internals + integer(I4B), dimension(:), allocatable :: ind + integer(I4B) :: direction + + if (self%nbody == 0) return - call swiftest_util_sort_qsort_DP(arr) + if (ascending) then + direction = 1 + else + direction = -1 + end if - return - end subroutine swiftest_util_sort_dp + associate(pl => self, npl => self%nbody) + select case(sortby) + case("Gmass","mass") + call util_sort(direction * pl%Gmass(1:npl), ind) + case("rhill") + call util_sort(direction * pl%rhill(1:npl), ind) + case("renc") + call util_sort(direction * pl%renc(1:npl), ind) + case("radius") + call util_sort(direction * pl%radius(1:npl), ind) + case("density") + call util_sort(direction * pl%density(1:npl), ind) + case("k2") + call util_sort(direction * pl%k2(1:npl), ind) + case("Q") + call util_sort(direction * pl%Q(1:npl), ind) + case("tlag") + call util_sort(direction * pl%tlag(1:npl), ind) + case("nplenc") + call util_sort(direction * pl%nplenc(1:npl), ind) + case("ntpenc") + call util_sort(direction * pl%ntpenc(1:npl), ind) + case("lmtiny", "nplm", "nplplm", "kin", "rbeg", "rend", "vbeg", "Ip", "rot", "k_plpl", "nplpl") + write(*,*) 'Cannot sort by ' // trim(adjustl(sortby)) // '. Component not sortable!' + case default ! Look for components in the parent class + call swiftest_util_sort_body(pl, sortby, ascending) + return + end select + call pl%rearrange(ind) - pure module subroutine swiftest_util_sort_index_dp(arr, ind) - !! author: David A. Minton - !! - !! Sort input DP precision array by index in ascending numerical order using quick sort. - !! This algorithm works well for partially sorted arrays (which is usually the case here). - !! If ind is supplied already allocated, we assume it is an existing index array (e.g. a previously - !! sorted array). If it is not allocated, this subroutine swiftest_allocates it. - !! - implicit none - ! Arguments - real(DP), dimension(:), intent(in) :: arr - integer(I4B), dimension(:), allocatable, intent(inout) :: ind - ! Internals - integer(I4B) :: n, i - real(DP), dimension(:), allocatable :: tmparr + end associate - n = size(arr) - if (.not.allocated(ind)) then - allocate(ind(n)) - ind = [(i, i=1, n)] - end if - allocate(tmparr, mold=arr) - tmparr(:) = arr(ind(:)) - call swiftest_util_sort_qsort_DP(tmparr, ind) - return - end subroutine swiftest_util_sort_index_dp + end subroutine swiftest_util_sort_pl - recursive pure subroutine swiftest_util_sort_qsort_DP(arr, ind) + module subroutine swiftest_util_sort_tp(self, sortby, ascending) !! author: David A. Minton !! - !! Sort input DP precision array by index in ascending numerical order using quicksort sort. - !! + !! Sort a Swiftest test particle object in-place. + !! sortby is a string indicating which array component to sort. implicit none ! Arguments - real(DP), dimension(:), intent(inout) :: arr - integer(I4B),dimension(:),intent(out), optional :: ind - !! Internals - integer :: iq - - if (size(arr) > 1) then - if (present(ind)) then - call swiftest_util_sort_partition_DP(arr, iq, ind) - call swiftest_util_sort_qsort_DP(arr(:iq-1),ind(:iq-1)) - call swiftest_util_sort_qsort_DP(arr(iq:), ind(iq:)) - else - call swiftest_util_sort_partition_DP(arr, iq) - call swiftest_util_sort_qsort_DP(arr(:iq-1)) - call swiftest_util_sort_qsort_DP(arr(iq:)) - end if + class(swiftest_tp), intent(inout) :: self !! Swiftest test particle object + character(*), intent(in) :: sortby !! Sorting attribute + logical, intent(in) :: ascending !! Logical flag indicating whether or not the sorting should be in ascending or descending order + ! Internals + integer(I4B), dimension(:), allocatable :: ind + integer(I4B) :: direction + + if (self%nbody == 0) return + + if (ascending) then + direction = 1 + else + direction = -1 end if - return - end subroutine swiftest_util_sort_qsort_DP + associate(tp => self, ntp => self%nbody) + select case(sortby) + case("nplenc") + call util_sort(direction * tp%nplenc(1:ntp), ind) + case default ! Look for components in the parent class + call swiftest_util_sort_body(tp, sortby, ascending) + return + end select - - pure subroutine swiftest_util_sort_partition_DP(arr, marker, ind) - !! author: David A. Minton - !! - !! Partition function for quicksort on DP type - !! - implicit none - ! Arguments - real(DP), intent(inout), dimension(:) :: arr - integer(I4B), intent(inout), dimension(:), optional :: ind - integer(I4B), intent(out) :: marker - ! Internals - integer(I4B) :: i, j, itmp, narr, ipiv - real(DP) :: temp - real(DP) :: x ! pivot point + call tp%rearrange(ind) - narr = size(arr) + end associate - ! Get center as pivot, as this is likely partially sorted - ipiv = narr / 2 - x = arr(ipiv) - i = 0 - j = narr + 1 - - do - j = j - 1 - do - if (arr(j) <= x) exit - j = j - 1 - end do - i = i + 1 - do - if (arr(i) >= x) exit - i = i + 1 - end do - if (i < j) then - ! exchange A(i) and A(j) - temp = arr(i) - arr(i) = arr(j) - arr(j) = temp - if (present(ind)) then - itmp = ind(i) - ind(i) = ind(j) - ind(j) = itmp - end if - else if (i == j) then - marker = i + 1 - return - else - marker = i - return - endif - end do - return - end subroutine swiftest_util_sort_partition_DP - + end subroutine swiftest_util_sort_tp - pure module subroutine swiftest_util_sort_i4b(arr) + module subroutine swiftest_util_sort_rearrange_body(self, ind) !! author: David A. Minton !! - !! Sort input integer array in place into ascending numerical order using quick sort. - !! This algorithm works well for partially sorted arrays (which is usually the case here) - !! + !! Rearrange Swiftest body structure in-place from an index list. + !! This is a helper utility used to make polymorphic sorting work on Swiftest structures. implicit none ! Arguments - integer(I4B), dimension(:), intent(inout) :: arr + class(swiftest_body), intent(inout) :: self !! Swiftest body object + integer(I4B), dimension(:), intent(in) :: ind !! Index array used to restructure the body (should contain all 1:n index values in the desired order) - call swiftest_util_sort_qsort_I4B(arr) + associate(n => self%nbody) + call util_sort_rearrange(self%id, ind, n) + call util_sort_rearrange(self%lmask, ind, n) + call util_sort_rearrange(self%info, ind, n) + call util_sort_rearrange(self%status, ind, n) + call util_sort_rearrange(self%ldiscard, ind, n) + call util_sort_rearrange(self%lcollision, ind, n) + call util_sort_rearrange(self%lencounter, ind, n) + call util_sort_rearrange(self%rh, ind, n) + call util_sort_rearrange(self%vh, ind, n) + call util_sort_rearrange(self%rb, ind, n) + call util_sort_rearrange(self%vb, ind, n) + call util_sort_rearrange(self%ah, ind, n) + call util_sort_rearrange(self%aobl, ind, n) + call util_sort_rearrange(self%agr, ind, n) + call util_sort_rearrange(self%atide, ind, n) + call util_sort_rearrange(self%ir3h, ind, n) + call util_sort_rearrange(self%isperi, ind, n) + call util_sort_rearrange(self%peri, ind, n) + call util_sort_rearrange(self%atp, ind, n) + call util_sort_rearrange(self%mu, ind, n) + call util_sort_rearrange(self%a, ind, n) + call util_sort_rearrange(self%e, ind, n) + call util_sort_rearrange(self%inc, ind, n) + call util_sort_rearrange(self%capom, ind, n) + call util_sort_rearrange(self%omega, ind, n) + call util_sort_rearrange(self%capm, ind, n) + end associate return - end subroutine swiftest_util_sort_i4b + end subroutine swiftest_util_sort_rearrange_body - pure module subroutine swiftest_util_sort_index_I4B(arr, ind) + module subroutine swiftest_util_sort_rearrange_arr_info(arr, ind, n) !! author: David A. Minton !! - !! Sort input integer array by index in ascending numerical order using quicksort. - !! If ind is supplied already allocated, we assume it is an existing index array (e.g. a previously - !! sorted array). If it is not allocated, this subroutine swiftest_allocates it. - !! + !! Rearrange a single array of particle information type in-place from an index list. implicit none ! Arguments - integer(I4B), dimension(:), intent(in) :: arr - integer(I4B), dimension(:), allocatable, intent(inout) :: ind + type(swiftest_particle_info), dimension(:), allocatable, intent(inout) :: arr !! Destination array + integer(I4B), dimension(:), intent(in) :: ind !! Index to rearrange against + integer(I4B), intent(in) :: n !! Number of elements in arr and ind to rearrange ! Internals - integer(I4B) :: n, i - integer(I4B), dimension(:), allocatable :: tmparr - - n = size(arr) - if (.not.allocated(ind)) then - allocate(ind(n)) - ind = [(i, i=1, n)] - end if - allocate(tmparr, mold=arr) - tmparr(:) = arr(ind(:)) - call swiftest_util_sort_qsort_I4B(tmparr, ind) + type(swiftest_particle_info), dimension(:), allocatable :: tmp !! Temporary copy of array used during rearrange operation - return - end subroutine swiftest_util_sort_index_I4B - - - pure module subroutine swiftest_util_sort_index_I4B_I8Bind(arr, ind) - !! author: David A. Minton - !! - !! Sort input integer array by index in ascending numerical order using quicksort. - !! If ind is supplied already allocated, we assume it is an existing index array (e.g. a previously - !! sorted array). If it is not allocated, this subroutine swiftest_allocates it. - !! - implicit none - ! Arguments - integer(I4B), dimension(:), intent(in) :: arr - integer(I8B), dimension(:), allocatable, intent(inout) :: ind - ! Internals - integer(I8B) :: n, i - integer(I4B), dimension(:), allocatable :: tmparr - - n = size(arr) - if (.not.allocated(ind)) then - allocate(ind(n)) - ind = [(i, i=1_I8B, n)] - end if - allocate(tmparr, mold=arr) - tmparr(:) = arr(ind(:)) - call swiftest_util_sort_qsort_I4B_I8Bind(tmparr, ind) - - return - end subroutine swiftest_util_sort_index_I4B_I8Bind - - - recursive pure subroutine swiftest_util_sort_qsort_I4B(arr, ind) - !! author: David A. Minton - !! - !! Sort input I4B array by index in ascending numerical order using quicksort. - !! - implicit none - ! Arguments - integer(I4B), dimension(:), intent(inout) :: arr - integer(I4B), dimension(:), intent(out), optional :: ind - ! Internals - integer(I4B) :: iq - - if (size(arr) > 1) then - if (present(ind)) then - call swiftest_util_sort_partition_I4B(arr, iq, ind) - call swiftest_util_sort_qsort_I4B(arr(:iq-1),ind(:iq-1)) - call swiftest_util_sort_qsort_I4B(arr(iq:), ind(iq:)) - else - call swiftest_util_sort_partition_I4B(arr, iq) - call swiftest_util_sort_qsort_I4B(arr(:iq-1)) - call swiftest_util_sort_qsort_I4B(arr(iq:)) - end if - end if - - return - end subroutine swiftest_util_sort_qsort_I4B - - - recursive pure subroutine swiftest_util_sort_qsort_I4B_I8Bind(arr, ind) - !! author: David A. Minton - !! - !! Sort input I4B array by index in ascending numerical order using quicksort. - !! - implicit none - ! Arguments - integer(I4B), dimension(:), intent(inout) :: arr - integer(I8B), dimension(:), intent(out), optional :: ind - ! Internals - integer(I8B) :: iq - - if (size(arr) > 1_I8B) then - if (present(ind)) then - call swiftest_util_sort_partition_I4B_I8Bind(arr, iq, ind) - call swiftest_util_sort_qsort_I4B_I8Bind(arr(:iq-1_I8B),ind(:iq-1_I8B)) - call swiftest_util_sort_qsort_I4B_I8Bind(arr(iq:), ind(iq:)) - else - call swiftest_util_sort_partition_I4B_I8Bind(arr, iq) - call swiftest_util_sort_qsort_I4B_I8Bind(arr(:iq-1_I8B)) - call swiftest_util_sort_qsort_I4B_I8Bind(arr(iq:)) - end if - end if - - return - end subroutine swiftest_util_sort_qsort_I4B_I8Bind - - - recursive pure subroutine swiftest_util_sort_qsort_I8B_I8Bind(arr, ind) - !! author: David A. Minton - !! - !! Sort input I8B array by index in ascending numerical order using quicksort. - !! - implicit none - ! Arguments - integer(I8B), dimension(:), intent(inout) :: arr - integer(I8B), dimension(:), intent(out), optional :: ind - ! Internals - integer(I8B) :: iq - - if (size(arr) > 1_I8B) then - if (present(ind)) then - call swiftest_util_sort_partition_I8B_I8Bind(arr, iq, ind) - call swiftest_util_sort_qsort_I8B_I8Bind(arr(:iq-1_I8B),ind(:iq-1_I8B)) - call swiftest_util_sort_qsort_I8B_I8Bind(arr(iq:), ind(iq:)) - else - call swiftest_util_sort_partition_I8B_I8Bind(arr, iq) - call swiftest_util_sort_qsort_I8B_I8Bind(arr(:iq-1_I8B)) - call swiftest_util_sort_qsort_I8B_I8Bind(arr(iq:)) - end if - end if - - return - end subroutine swiftest_util_sort_qsort_I8B_I8Bind - - - pure subroutine swiftest_util_sort_partition_I4B(arr, marker, ind) - !! author: David A. Minton - !! - !! Partition function for quicksort on I4B type - !! - implicit none - ! Arguments - integer(I4B), intent(inout), dimension(:) :: arr - integer(I4B), intent(inout), dimension(:), optional :: ind - integer(I4B), intent(out) :: marker - ! Internals - integer(I4B) :: i, j, itmp, narr, ipiv - integer(I4B) :: temp - integer(I4B) :: x ! pivot point - - narr = size(arr) - - ! Get center as pivot, as this is likely partially sorted - ipiv = narr / 2 - x = arr(ipiv) - i = 0 - j = narr + 1 - - do - j = j - 1 - do - if (arr(j) <= x) exit - j = j - 1 - end do - i = i + 1 - do - if (arr(i) >= x) exit - i = i + 1 - end do - if (i < j) then - ! exchange A(i) and A(j) - temp = arr(i) - arr(i) = arr(j) - arr(j) = temp - if (present(ind)) then - itmp = ind(i) - ind(i) = ind(j) - ind(j) = itmp - end if - else if (i == j) then - marker = i + 1 - return - else - marker = i - return - endif - end do - - return - end subroutine swiftest_util_sort_partition_I4B - - - pure subroutine swiftest_util_sort_partition_I4B_I8Bind(arr, marker, ind) - !! author: David A. Minton - !! - !! Partition function for quicksort on I4B type - !! - implicit none - ! Arguments - integer(I4B), intent(inout), dimension(:) :: arr - integer(I8B), intent(inout), dimension(:), optional :: ind - integer(I8B), intent(out) :: marker - ! Internals - integer(I8B) :: i, j, itmp, narr, ipiv - integer(I4B) :: temp - integer(I8B) :: x ! pivot point - - narr = size(arr) - - ! Get center as pivot, as this is likely partially sorted - ipiv = narr / 2_I8B - x = arr(ipiv) - i = 0_I8B - j = narr + 1_I8B - - do - j = j - 1_I8B - do - if (arr(j) <= x) exit - j = j - 1_I8B - end do - i = i + 1_I8B - do - if (arr(i) >= x) exit - i = i + 1_I8B - end do - if (i < j) then - ! exchange A(i) and A(j) - temp = arr(i) - arr(i) = arr(j) - arr(j) = temp - if (present(ind)) then - itmp = ind(i) - ind(i) = ind(j) - ind(j) = itmp - end if - else if (i == j) then - marker = i + 1_I8B - return - else - marker = i - return - endif - end do - - return - end subroutine swiftest_util_sort_partition_I4B_I8Bind - - - pure subroutine swiftest_util_sort_partition_I8B_I8Bind(arr, marker, ind) - !! author: David A. Minton - !! - !! Partition function for quicksort on I8B type with I8B index - !! - implicit none - ! Arguments - integer(I8B), intent(inout), dimension(:) :: arr - integer(I8B), intent(inout), dimension(:), optional :: ind - integer(I8B), intent(out) :: marker - ! Internals - integer(I8B) :: i, j, itmp, narr, ipiv - integer(I8B) :: temp - integer(I8B) :: x ! pivot point - - narr = size(arr) - - ! Get center as pivot, as this is likely partially sorted - ipiv = narr / 2_I8B - x = arr(ipiv) - i = 0_I8B - j = narr + 1_I8B - - do - j = j - 1_I8B - do - if (arr(j) <= x) exit - j = j - 1_I8B - end do - i = i + 1_I8B - do - if (arr(i) >= x) exit - i = i + 1_I8B - end do - if (i < j) then - ! exchange A(i) and A(j) - temp = arr(i) - arr(i) = arr(j) - arr(j) = temp - if (present(ind)) then - itmp = ind(i) - ind(i) = ind(j) - ind(j) = itmp - end if - else if (i == j) then - marker = i + 1_I8B - return - else - marker = i - return - endif - end do - - return - end subroutine swiftest_util_sort_partition_I8B_I8Bind - - - pure module subroutine swiftest_util_sort_sp(arr) - !! author: David A. Minton - !! - !! Sort input DP precision array in place into ascending numerical order using quicksort. - !! - implicit none - ! Arguments - real(SP), dimension(:), intent(inout) :: arr - - call swiftest_util_sort_qsort_SP(arr) - - return - end subroutine swiftest_util_sort_sp - - - pure module subroutine swiftest_util_sort_index_sp(arr, ind) - !! author: David A. Minton - !! - !! Sort input DP precision array by index in ascending numerical order using quicksort. - !! If ind is supplied already allocated, we assume it is an existing index array (e.g. a previously - !! sorted array). If it is not allocated, this subroutine swiftest_allocates it. - !! - implicit none - ! Arguments - real(SP), dimension(:), intent(in) :: arr - integer(I4B), dimension(:), allocatable, intent(inout) :: ind - ! Internals - integer(I4B) :: n, i - real(SP), dimension(:), allocatable :: tmparr - - n = size(arr) - if (.not.allocated(ind)) then - allocate(ind(n)) - ind = [(i, i=1, n)] - end if - allocate(tmparr, mold=arr) - tmparr(:) = arr(ind(:)) - call swiftest_util_sort_qsort_SP(tmparr, ind) - - return - end subroutine swiftest_util_sort_index_sp - - - recursive pure subroutine swiftest_util_sort_qsort_SP(arr, ind) - !! author: David A. Minton - !! - !! Sort input DP precision array by index in ascending numerical order using quicksort. - !! - implicit none - ! Arguments - real(SP), dimension(:), intent(inout) :: arr - integer(I4B),dimension(:),intent(out), optional :: ind - !! Internals - integer :: iq - - if (size(arr) > 1) then - if (present(ind)) then - call swiftest_util_sort_partition_SP(arr, iq, ind) - call swiftest_util_sort_qsort_SP(arr(:iq-1),ind(:iq-1)) - call swiftest_util_sort_qsort_SP(arr(iq:), ind(iq:)) - else - call swiftest_util_sort_partition_SP(arr, iq) - call swiftest_util_sort_qsort_SP(arr(:iq-1)) - call swiftest_util_sort_qsort_SP(arr(iq:)) - end if - end if - - return - end subroutine swiftest_util_sort_qsort_SP - - - pure subroutine swiftest_util_sort_partition_SP(arr, marker, ind) - !! author: David A. Minton - !! - !! Partition function for quicksort on SP type - !! - implicit none - ! Arguments - real(SP), intent(inout), dimension(:) :: arr - integer(I4B), intent(inout), dimension(:), optional :: ind - integer(I4B), intent(out) :: marker - ! Internals - integer(I4B) :: i, j, itmp, narr, ipiv - real(SP) :: temp - real(SP) :: x ! pivot point - - narr = size(arr) - - ! Get center as pivot, as this is likely partially sorted - ipiv = narr / 2 - x = arr(ipiv) - i = 0 - j = narr + 1 - - do - j = j - 1 - do - if (arr(j) <= x) exit - j = j - 1 - end do - i = i + 1 - do - if (arr(i) >= x) exit - i = i + 1 - end do - if (i < j) then - ! exchange A(i) and A(j) - temp = arr(i) - arr(i) = arr(j) - arr(j) = temp - if (present(ind)) then - itmp = ind(i) - ind(i) = ind(j) - ind(j) = itmp - end if - else if (i == j) then - marker = i + 1 - return - else - marker = i - return - endif - end do - - return - end subroutine swiftest_util_sort_partition_SP - - - module subroutine swiftest_util_sort_pl(self, sortby, ascending) - !! author: David A. Minton - !! - !! Sort a Swiftest massive body object in-place. - !! sortby is a string indicating which array component to sort. - implicit none - ! Arguments - class(swiftest_pl), intent(inout) :: self !! Swiftest massive body object - character(*), intent(in) :: sortby !! Sorting attribute - logical, intent(in) :: ascending !! Logical flag indicating whether or not the sorting should be in ascending or descending order - ! Internals - integer(I4B), dimension(:), allocatable :: ind - integer(I4B) :: direction - - if (self%nbody == 0) return - - if (ascending) then - direction = 1 - else - direction = -1 - end if - - associate(pl => self, npl => self%nbody) - select case(sortby) - case("Gmass","mass") - call swiftest_util_sort(direction * pl%Gmass(1:npl), ind) - case("rhill") - call swiftest_util_sort(direction * pl%rhill(1:npl), ind) - case("renc") - call swiftest_util_sort(direction * pl%renc(1:npl), ind) - case("radius") - call swiftest_util_sort(direction * pl%radius(1:npl), ind) - case("density") - call swiftest_util_sort(direction * pl%density(1:npl), ind) - case("k2") - call swiftest_util_sort(direction * pl%k2(1:npl), ind) - case("Q") - call swiftest_util_sort(direction * pl%Q(1:npl), ind) - case("tlag") - call swiftest_util_sort(direction * pl%tlag(1:npl), ind) - case("nplenc") - call swiftest_util_sort(direction * pl%nplenc(1:npl), ind) - case("ntpenc") - call swiftest_util_sort(direction * pl%ntpenc(1:npl), ind) - case("lmtiny", "nplm", "nplplm", "kin", "rbeg", "rend", "vbeg", "Ip", "rot", "k_plpl", "nplpl") - write(*,*) 'Cannot sort by ' // trim(adjustl(sortby)) // '. Component not sortable!' - case default ! Look for components in the parent class - call swiftest_util_sort_body(pl, sortby, ascending) - return - end select - - call pl%rearrange(ind) - - end associate - - return - end subroutine swiftest_util_sort_pl - - - module subroutine swiftest_util_sort_tp(self, sortby, ascending) - !! author: David A. Minton - !! - !! Sort a Swiftest test particle object in-place. - !! sortby is a string indicating which array component to sort. - implicit none - ! Arguments - class(swiftest_tp), intent(inout) :: self !! Swiftest test particle object - character(*), intent(in) :: sortby !! Sorting attribute - logical, intent(in) :: ascending !! Logical flag indicating whether or not the sorting should be in ascending or descending order - ! Internals - integer(I4B), dimension(:), allocatable :: ind - integer(I4B) :: direction - - if (self%nbody == 0) return - - if (ascending) then - direction = 1 - else - direction = -1 - end if - - associate(tp => self, ntp => self%nbody) - select case(sortby) - case("nplenc") - call swiftest_util_sort(direction * tp%nplenc(1:ntp), ind) - case default ! Look for components in the parent class - call swiftest_util_sort_body(tp, sortby, ascending) - return - end select - - call tp%rearrange(ind) - - end associate - - return - end subroutine swiftest_util_sort_tp - - - module subroutine swiftest_util_sort_rearrange_body(self, ind) - !! author: David A. Minton - !! - !! Rearrange Swiftest body structure in-place from an index list. - !! This is a helper utility used to make polymorphic sorting work on Swiftest structures. - implicit none - ! Arguments - class(swiftest_body), intent(inout) :: self !! Swiftest body object - integer(I4B), dimension(:), intent(in) :: ind !! Index array used to restructure the body (should contain all 1:n index values in the desired order) - - associate(n => self%nbody) - call swiftest_util_sort_rearrange(self%id, ind, n) - call swiftest_util_sort_rearrange(self%lmask, ind, n) - call swiftest_util_sort_rearrange(self%info, ind, n) - call swiftest_util_sort_rearrange(self%status, ind, n) - call swiftest_util_sort_rearrange(self%ldiscard, ind, n) - call swiftest_util_sort_rearrange(self%lcollision, ind, n) - call swiftest_util_sort_rearrange(self%lencounter, ind, n) - call swiftest_util_sort_rearrange(self%rh, ind, n) - call swiftest_util_sort_rearrange(self%vh, ind, n) - call swiftest_util_sort_rearrange(self%rb, ind, n) - call swiftest_util_sort_rearrange(self%vb, ind, n) - call swiftest_util_sort_rearrange(self%ah, ind, n) - call swiftest_util_sort_rearrange(self%aobl, ind, n) - call swiftest_util_sort_rearrange(self%agr, ind, n) - call swiftest_util_sort_rearrange(self%atide, ind, n) - call swiftest_util_sort_rearrange(self%ir3h, ind, n) - call swiftest_util_sort_rearrange(self%isperi, ind, n) - call swiftest_util_sort_rearrange(self%peri, ind, n) - call swiftest_util_sort_rearrange(self%atp, ind, n) - call swiftest_util_sort_rearrange(self%mu, ind, n) - call swiftest_util_sort_rearrange(self%a, ind, n) - call swiftest_util_sort_rearrange(self%e, ind, n) - call swiftest_util_sort_rearrange(self%inc, ind, n) - call swiftest_util_sort_rearrange(self%capom, ind, n) - call swiftest_util_sort_rearrange(self%omega, ind, n) - call swiftest_util_sort_rearrange(self%capm, ind, n) - end associate - - return - end subroutine swiftest_util_sort_rearrange_body - - - pure module subroutine swiftest_util_sort_rearrange_arr_char_string(arr, ind, n) - !! author: David A. Minton - !! - !! Rearrange a single array of character string in-place from an index list. - implicit none - ! Arguments - character(len=STRMAX), dimension(:), allocatable, intent(inout) :: arr !! Destination array - integer(I4B), dimension(:), intent(in) :: ind !! Index to rearrange against - integer(I4B), intent(in) :: n !! Number of elements in arr and ind to rearrange - ! Internals - character(len=STRMAX), dimension(:), allocatable :: tmp !! Temporary copy of arry used during rearrange operation - - if (.not. allocated(arr) .or. n <= 0) return - allocate(tmp, mold=arr) - tmp(1:n) = arr(ind) - call move_alloc(tmp, arr) - - return - end subroutine swiftest_util_sort_rearrange_arr_char_string - - - pure module subroutine swiftest_util_sort_rearrange_arr_DP(arr, ind, n) - !! author: David A. Minton - !! - !! Rearrange a single array of DP type in-place from an index list. - implicit none - ! Arguments - real(DP), dimension(:), allocatable, intent(inout) :: arr !! Destination array - integer(I4B), dimension(:), intent(in) :: ind !! Index to rearrange against - integer(I4B), intent(in) :: n !! Number of elements in arr and ind to rearrange - ! Internals - real(DP), dimension(:), allocatable :: tmp !! Temporary copy of array used during rearrange operation - - if (.not. allocated(arr) .or. n <= 0) return - allocate(tmp, mold=arr) - tmp(1:n) = arr(ind) - call move_alloc(tmp, arr) - - return - end subroutine swiftest_util_sort_rearrange_arr_DP - - - pure module subroutine swiftest_util_sort_rearrange_arr_DPvec(arr, ind, n) - !! author: David A. Minton - !! - !! Rearrange a single array of (NDIM,n) DP-type vectors in-place from an index list. - implicit none - ! Arguments - real(DP), dimension(:,:), allocatable, intent(inout) :: arr !! Destination array - integer(I4B), dimension(:), intent(in) :: ind !! Index to rearrange against - integer(I4B), intent(in) :: n !! Number of elements in arr and ind to rearrange - ! Internals - real(DP), dimension(:,:), allocatable :: tmp !! Temporary copy of array used during rearrange operation - - if (.not. allocated(arr) .or. n <= 0) return - allocate(tmp, mold=arr) - tmp(:,1:n) = arr(:, ind) - call move_alloc(tmp, arr) - - return - end subroutine swiftest_util_sort_rearrange_arr_DPvec - - - pure module subroutine swiftest_util_sort_rearrange_arr_I4B(arr, ind, n) - !! author: David A. Minton - !! - !! Rearrange a single array of integers in-place from an index list. - implicit none - ! Arguments - integer(I4B), dimension(:), allocatable, intent(inout) :: arr !! Destination array - integer(I4B), dimension(:), intent(in) :: ind !! Index to rearrange against - integer(I4B), intent(in) :: n !! Number of elements in arr and ind to rearrange - ! Internals - integer(I4B), dimension(:), allocatable :: tmp !! Temporary copy of array used during rearrange operation - - if (.not. allocated(arr) .or. n <= 0) return - allocate(tmp, mold=arr) - tmp(1:n) = arr(ind) - call move_alloc(tmp, arr) - - return - end subroutine swiftest_util_sort_rearrange_arr_I4B - - pure module subroutine swiftest_util_sort_rearrange_arr_I4B_I8Bind(arr, ind, n) - !! author: David A. Minton - !! - !! Rearrange a single array of integers in-place from an index list. - implicit none - ! Arguments - integer(I4B), dimension(:), allocatable, intent(inout) :: arr !! Destination array - integer(I8B), dimension(:), intent(in) :: ind !! Index to rearrange against - integer(I8B), intent(in) :: n !! Number of elements in arr and ind to rearrange - ! Internals - integer(I4B), dimension(:), allocatable :: tmp !! Temporary copy of array used during rearrange operation - - if (.not. allocated(arr) .or. n <= 0_I8B) return - allocate(tmp, mold=arr) - tmp(1:n) = arr(ind) - call move_alloc(tmp, arr) - - return - end subroutine swiftest_util_sort_rearrange_arr_I4B_I8Bind - - - module subroutine swiftest_util_sort_rearrange_arr_info(arr, ind, n) - !! author: David A. Minton - !! - !! Rearrange a single array of particle information type in-place from an index list. - implicit none - ! Arguments - type(swiftest_particle_info), dimension(:), allocatable, intent(inout) :: arr !! Destination array - integer(I4B), dimension(:), intent(in) :: ind !! Index to rearrange against - integer(I4B), intent(in) :: n !! Number of elements in arr and ind to rearrange - ! Internals - type(swiftest_particle_info), dimension(:), allocatable :: tmp !! Temporary copy of array used during rearrange operation - - if (.not. allocated(arr) .or. n <= 0) return - allocate(tmp, mold=arr) + if (.not. allocated(arr) .or. n <= 0) return + allocate(tmp, mold=arr) call swiftest_util_copy_particle_info_arr(arr, tmp, ind) call move_alloc(tmp, arr) return end subroutine swiftest_util_sort_rearrange_arr_info - + pure module subroutine swiftest_util_sort_rearrange_arr_kin(arr, ind, n) !! author: David A. Minton @@ -4067,48 +2911,6 @@ pure module subroutine swiftest_util_sort_rearrange_arr_kin(arr, ind, n) end subroutine swiftest_util_sort_rearrange_arr_kin - pure module subroutine swiftest_util_sort_rearrange_arr_logical(arr, ind, n) - !! author: David A. Minton - !! - !! Rearrange a single array of logicals in-place from an index list. - implicit none - ! Arguments - logical, dimension(:), allocatable, intent(inout) :: arr !! Destination array - integer(I4B), dimension(:), intent(in) :: ind !! Index to rearrange against - integer(I4B), intent(in) :: n !! Number of elements in arr and ind to rearrange - ! Internals - logical, dimension(:), allocatable :: tmp !! Temporary copy of array used during rearrange operation - - if (.not. allocated(arr) .or. n <= 0) return - allocate(tmp, mold=arr) - tmp(1:n) = arr(ind) - call move_alloc(tmp, arr) - - return - end subroutine swiftest_util_sort_rearrange_arr_logical - - - pure module subroutine swiftest_util_sort_rearrange_arr_logical_I8Bind(arr, ind, n) - !! author: David A. Minton - !! - !! Rearrange a single array of logicals in-place from an index list. - implicit none - ! Arguments - logical, dimension(:), allocatable, intent(inout) :: arr !! Destination array - integer(I8B), dimension(:), intent(in) :: ind !! Index to rearrange against - integer(I8B), intent(in) :: n !! Number of elements in arr and ind to rearrange - ! Internals - logical, dimension(:), allocatable :: tmp !! Temporary copy of array used during rearrange operation - - if (.not. allocated(arr) .or. n <= 0) return - allocate(tmp, mold=arr) - tmp(1:n) = arr(ind) - call move_alloc(tmp, arr) - - return - end subroutine swiftest_util_sort_rearrange_arr_logical_I8Bind - - module subroutine swiftest_util_sort_rearrange_pl(self, ind) !! author: David A. Minton !! @@ -4119,23 +2921,23 @@ module subroutine swiftest_util_sort_rearrange_pl(self, ind) integer(I4B), dimension(:), intent(in) :: ind !! Index array used to restructure the body (should contain all 1:n index values in the desired order) associate(pl => self, npl => self%nbody) - call swiftest_util_sort_rearrange(pl%mass, ind, npl) - call swiftest_util_sort_rearrange(pl%Gmass, ind, npl) - call swiftest_util_sort_rearrange(pl%rhill, ind, npl) - call swiftest_util_sort_rearrange(pl%renc, ind, npl) - call swiftest_util_sort_rearrange(pl%radius, ind, npl) - call swiftest_util_sort_rearrange(pl%density, ind, npl) - call swiftest_util_sort_rearrange(pl%rbeg, ind, npl) - call swiftest_util_sort_rearrange(pl%vbeg, ind, npl) - call swiftest_util_sort_rearrange(pl%Ip, ind, npl) - call swiftest_util_sort_rearrange(pl%rot, ind, npl) - call swiftest_util_sort_rearrange(pl%k2, ind, npl) - call swiftest_util_sort_rearrange(pl%Q, ind, npl) - call swiftest_util_sort_rearrange(pl%tlag, ind, npl) - call swiftest_util_sort_rearrange(pl%kin, ind, npl) - call swiftest_util_sort_rearrange(pl%lmtiny, ind, npl) - call swiftest_util_sort_rearrange(pl%nplenc, ind, npl) - call swiftest_util_sort_rearrange(pl%ntpenc, ind, npl) + call util_sort_rearrange(pl%mass, ind, npl) + call util_sort_rearrange(pl%Gmass, ind, npl) + call util_sort_rearrange(pl%rhill, ind, npl) + call util_sort_rearrange(pl%renc, ind, npl) + call util_sort_rearrange(pl%radius, ind, npl) + call util_sort_rearrange(pl%density, ind, npl) + call util_sort_rearrange(pl%rbeg, ind, npl) + call util_sort_rearrange(pl%vbeg, ind, npl) + call util_sort_rearrange(pl%Ip, ind, npl) + call util_sort_rearrange(pl%rot, ind, npl) + call util_sort_rearrange(pl%k2, ind, npl) + call util_sort_rearrange(pl%Q, ind, npl) + call util_sort_rearrange(pl%tlag, ind, npl) + call util_sort_rearrange(pl%kin, ind, npl) + call util_sort_rearrange(pl%lmtiny, ind, npl) + call util_sort_rearrange(pl%nplenc, ind, npl) + call util_sort_rearrange(pl%ntpenc, ind, npl) if (allocated(pl%k_plpl)) deallocate(pl%k_plpl) @@ -4157,7 +2959,7 @@ module subroutine swiftest_util_sort_rearrange_tp(self, ind) integer(I4B), dimension(:), intent(in) :: ind !! Index array used to restructure the body (should contain all 1:n index values in the desired order) associate(tp => self, ntp => self%nbody) - call swiftest_util_sort_rearrange(tp%nplenc, ind, ntp) + call util_sort_rearrange(tp%nplenc, ind, ntp) if (allocated(tp%k_pltp)) deallocate(tp%k_pltp) @@ -4168,220 +2970,6 @@ module subroutine swiftest_util_sort_rearrange_tp(self, ind) end subroutine swiftest_util_sort_rearrange_tp - module subroutine swiftest_util_spill_arr_char_string(keeps, discards, lspill_list, ldestructive) - !! author: David A. Minton - !! - !! Performs a spill operation on a single array of type character strings - !! This is the inverse of a spill operation - implicit none - ! Arguments - character(len=STRMAX), dimension(:), allocatable, intent(inout) :: keeps !! Array of values to keep - character(len=STRMAX), dimension(:), allocatable, intent(inout) :: discards !! Array of discards - logical, dimension(:), intent(in) :: lspill_list !! Logical array of bodies to spill into the discards - logical, intent(in) :: ldestructive !! Logical flag indicating whether or not this operation should alter the keeps array or not - ! Internals - integer(I4B) :: nspill, nkeep, nlist - character(len=STRMAX), dimension(:), allocatable :: tmp !! Array of values to keep - - nkeep = count(.not.lspill_list(:)) - nspill = count(lspill_list(:)) - nlist = size(lspill_list(:)) - - if (.not.allocated(keeps) .or. nspill == 0) return - if (.not.allocated(discards)) then - allocate(discards(nspill)) - else if (size(discards) /= nspill) then - deallocate(discards) - allocate(discards(nspill)) - end if - - discards(:) = pack(keeps(1:nlist), lspill_list(1:nlist)) - if (ldestructive) then - if (nkeep > 0) then - allocate(tmp(nkeep)) - tmp(:) = pack(keeps(1:nlist), .not. lspill_list(1:nlist)) - call move_alloc(tmp, keeps) - else - deallocate(keeps) - end if - end if - - return - end subroutine swiftest_util_spill_arr_char_string - - - module subroutine swiftest_util_spill_arr_DP(keeps, discards, lspill_list, ldestructive) - !! author: David A. Minton - !! - !! Performs a spill operation on a single array of type DP - !! This is the inverse of a spill operation - implicit none - ! Arguments - real(DP), dimension(:), allocatable, intent(inout) :: keeps !! Array of values to keep - real(DP), dimension(:), allocatable, intent(inout) :: discards !! Array of discards - logical, dimension(:), intent(in) :: lspill_list !! Logical array of bodies to spill into the discardss - logical, intent(in) :: ldestructive !! Logical flag indicating whether or not this operation should alter the keeps array or not - ! Internals - integer(I4B) :: nspill, nkeep, nlist - real(DP), dimension(:), allocatable :: tmp !! Array of values to keep - - nkeep = count(.not.lspill_list(:)) - nspill = count(lspill_list(:)) - nlist = size(lspill_list(:)) - - if (.not.allocated(keeps) .or. nspill == 0) return - if (.not.allocated(discards)) then - allocate(discards(nspill)) - else if (size(discards) /= nspill) then - deallocate(discards) - allocate(discards(nspill)) - end if - - discards(:) = pack(keeps(1:nlist), lspill_list(1:nlist)) - if (ldestructive) then - if (nkeep > 0) then - allocate(tmp(nkeep)) - tmp(:) = pack(keeps(1:nlist), .not. lspill_list(1:nlist)) - call move_alloc(tmp, keeps) - else - deallocate(keeps) - end if - end if - - return - end subroutine swiftest_util_spill_arr_DP - - - module subroutine swiftest_util_spill_arr_DPvec(keeps, discards, lspill_list, ldestructive) - !! author: David A. Minton - !! - !! Performs a spill operation on a single array of DP vectors with shape (NDIM, n) - !! This is the inverse of a spill operation - implicit none - ! Arguments - real(DP), dimension(:,:), allocatable, intent(inout) :: keeps !! Array of values to keep - real(DP), dimension(:,:), allocatable, intent(inout) :: discards !! Array discards - logical, dimension(:), intent(in) :: lspill_list !! Logical array of bodies to spill into the discards - logical, intent(in) :: ldestructive !! Logical flag indicating whether or not this operation should alter the keeps array or not - ! Internals - integer(I4B) :: i, nspill, nkeep, nlist - real(DP), dimension(:,:), allocatable :: tmp !! Array of values to keep - - nkeep = count(.not.lspill_list(:)) - nspill = count(lspill_list(:)) - nlist = size(lspill_list(:)) - - if (.not.allocated(keeps) .or. nspill == 0) return - if (.not.allocated(discards)) then - allocate(discards(NDIM, nspill)) - else if (size(discards, dim=2) /= nspill) then - deallocate(discards) - allocate(discards(NDIM, nspill)) - end if - - do i = 1, NDIM - discards(i,:) = pack(keeps(i,1:nlist), lspill_list(1:nlist)) - end do - if (ldestructive) then - if (nkeep > 0) then - allocate(tmp(NDIM, nkeep)) - do i = 1, NDIM - tmp(i, :) = pack(keeps(i, 1:nlist), .not. lspill_list(1:nlist)) - end do - call move_alloc(tmp, keeps) - else - deallocate(keeps) - end if - end if - - return - end subroutine swiftest_util_spill_arr_DPvec - - - module subroutine swiftest_util_spill_arr_I4B(keeps, discards, lspill_list, ldestructive) - !! author: David A. Minton - !! - !! Performs a spill operation on a single array of type I4B - !! This is the inverse of a spill operation - implicit none - ! Arguments - integer(I4B), dimension(:), allocatable, intent(inout) :: keeps !! Array of values to keep - integer(I4B), dimension(:), allocatable, intent(inout) :: discards !! Array of discards - logical, dimension(:), intent(in) :: lspill_list !! Logical array of bodies to spill into the discards - logical, intent(in) :: ldestructive !! Logical flag indicating whether or not this operation should alter the keeps array or not - ! Internals - integer(I4B) :: nspill, nkeep, nlist - integer(I4B), dimension(:), allocatable :: tmp !! Array of values to keep - - nkeep = count(.not.lspill_list(:)) - nspill = count(lspill_list(:)) - nlist = size(lspill_list(:)) - - if (.not.allocated(keeps) .or. nspill == 0) return - if (.not.allocated(discards)) then - allocate(discards(nspill)) - else if (size(discards) /= nspill) then - deallocate(discards) - allocate(discards(nspill)) - end if - - discards(:) = pack(keeps(1:nlist), lspill_list(1:nlist)) - if (ldestructive) then - if (nkeep > 0) then - allocate(tmp(nkeep)) - tmp(:) = pack(keeps(1:nlist), .not. lspill_list(1:nlist)) - call move_alloc(tmp, keeps) - else - deallocate(keeps) - end if - end if - - return - end subroutine swiftest_util_spill_arr_I4B - - - module subroutine swiftest_util_spill_arr_I8B(keeps, discards, lspill_list, ldestructive) - !! author: David A. Minton - !! - !! Performs a spill operation on a single array of type I4B - !! This is the inverse of a spill operation - implicit none - ! Arguments - integer(I8B), dimension(:), allocatable, intent(inout) :: keeps !! Array of values to keep - integer(I8B), dimension(:), allocatable, intent(inout) :: discards !! Array of discards - logical, dimension(:), intent(in) :: lspill_list !! Logical array of bodies to spill into the discards - logical, intent(in) :: ldestructive !! Logical flag indicating whether or not this operation should alter the keeps array or not - ! Internals - integer(I4B) :: nspill, nkeep, nlist - integer(I8B), dimension(:), allocatable :: tmp !! Array of values to keep - - nkeep = count(.not.lspill_list(:)) - nspill = count(lspill_list(:)) - nlist = size(lspill_list(:)) - - if (.not.allocated(keeps) .or. nspill == 0) return - if (.not.allocated(discards)) then - allocate(discards(nspill)) - else if (size(discards) /= nspill) then - deallocate(discards) - allocate(discards(nspill)) - end if - - discards(:) = pack(keeps(1:nlist), lspill_list(1:nlist)) - if (ldestructive) then - if (nkeep > 0) then - allocate(tmp(nkeep)) - tmp(:) = pack(keeps(1:nlist), .not. lspill_list(1:nlist)) - call move_alloc(tmp, keeps) - else - deallocate(keeps) - end if - end if - - return - end subroutine swiftest_util_spill_arr_I8B - - module subroutine swiftest_util_spill_arr_info(keeps, discards, lspill_list, ldestructive) !! author: David A. Minton !! @@ -4472,48 +3060,6 @@ module subroutine swiftest_util_spill_arr_kin(keeps, discards, lspill_list, ldes end subroutine swiftest_util_spill_arr_kin - module subroutine swiftest_util_spill_arr_logical(keeps, discards, lspill_list, ldestructive) - !! author: David A. Minton - !! - !! Performs a spill operation on a single array of logicals - !! This is the inverse of a spill operation - implicit none - ! Arguments - logical, dimension(:), allocatable, intent(inout) :: keeps !! Array of values to keep - logical, dimension(:), allocatable, intent(inout) :: discards !! Array of discards - logical, dimension(:), intent(in) :: lspill_list !! Logical array of bodies to spill into the discards - logical, intent(in) :: ldestructive !! Logical flag indicating whether or not this operation should alter the keeps array or no - ! Internals - integer(I4B) :: nspill, nkeep, nlist - logical, dimension(:), allocatable :: tmp !! Array of values to keep - - nkeep = count(.not.lspill_list(:)) - nspill = count(lspill_list(:)) - nlist = size(lspill_list(:)) - - if (.not.allocated(keeps) .or. nspill == 0) return - if (.not.allocated(discards)) then - allocate(discards(nspill)) - else if (size(discards) /= nspill) then - deallocate(discards) - allocate(discards(nspill)) - end if - - discards(:) = pack(keeps(1:nlist), lspill_list(1:nlist)) - if (ldestructive) then - if (nkeep > 0) then - allocate(tmp(nkeep)) - tmp(:) = pack(keeps(1:nlist), .not. lspill_list(1:nlist)) - call move_alloc(tmp, keeps) - else - deallocate(keeps) - end if - end if - - return - end subroutine swiftest_util_spill_arr_logical - - module subroutine swiftest_util_spill_body(self, discards, lspill_list, ldestructive) !! author: David A. Minton !! @@ -4532,32 +3078,32 @@ module subroutine swiftest_util_spill_body(self, discards, lspill_list, ldestruc !> Spill all the common components associate(keeps => self) - call swiftest_util_spill(keeps%id, discards%id, lspill_list, ldestructive) - call swiftest_util_spill(keeps%info, discards%info, lspill_list, ldestructive) - call swiftest_util_spill(keeps%lmask, discards%lmask, lspill_list, ldestructive) - call swiftest_util_spill(keeps%status, discards%status, lspill_list, ldestructive) - call swiftest_util_spill(keeps%ldiscard, discards%ldiscard, lspill_list, ldestructive) - call swiftest_util_spill(keeps%lcollision, discards%lcollision, lspill_list, ldestructive) - call swiftest_util_spill(keeps%lencounter, discards%lencounter, lspill_list, ldestructive) - call swiftest_util_spill(keeps%mu, discards%mu, lspill_list, ldestructive) - call swiftest_util_spill(keeps%rh, discards%rh, lspill_list, ldestructive) - call swiftest_util_spill(keeps%vh, discards%vh, lspill_list, ldestructive) - call swiftest_util_spill(keeps%rb, discards%rb, lspill_list, ldestructive) - call swiftest_util_spill(keeps%vb, discards%vb, lspill_list, ldestructive) - call swiftest_util_spill(keeps%ah, discards%ah, lspill_list, ldestructive) - call swiftest_util_spill(keeps%aobl, discards%aobl, lspill_list, ldestructive) - call swiftest_util_spill(keeps%agr, discards%agr, lspill_list, ldestructive) - call swiftest_util_spill(keeps%atide, discards%atide, lspill_list, ldestructive) - call swiftest_util_spill(keeps%ir3h, discards%ir3h, lspill_list, ldestructive) - call swiftest_util_spill(keeps%isperi, discards%isperi, lspill_list, ldestructive) - call swiftest_util_spill(keeps%peri, discards%peri, lspill_list, ldestructive) - call swiftest_util_spill(keeps%atp, discards%atp, lspill_list, ldestructive) - call swiftest_util_spill(keeps%a, discards%a, lspill_list, ldestructive) - call swiftest_util_spill(keeps%e, discards%e, lspill_list, ldestructive) - call swiftest_util_spill(keeps%inc, discards%inc, lspill_list, ldestructive) - call swiftest_util_spill(keeps%capom, discards%capom, lspill_list, ldestructive) - call swiftest_util_spill(keeps%omega, discards%omega, lspill_list, ldestructive) - call swiftest_util_spill(keeps%capm, discards%capm, lspill_list, ldestructive) + call util_spill(keeps%id, discards%id, lspill_list, ldestructive) + call util_spill(keeps%info, discards%info, lspill_list, ldestructive) + call util_spill(keeps%lmask, discards%lmask, lspill_list, ldestructive) + call util_spill(keeps%status, discards%status, lspill_list, ldestructive) + call util_spill(keeps%ldiscard, discards%ldiscard, lspill_list, ldestructive) + call util_spill(keeps%lcollision, discards%lcollision, lspill_list, ldestructive) + call util_spill(keeps%lencounter, discards%lencounter, lspill_list, ldestructive) + call util_spill(keeps%mu, discards%mu, lspill_list, ldestructive) + call util_spill(keeps%rh, discards%rh, lspill_list, ldestructive) + call util_spill(keeps%vh, discards%vh, lspill_list, ldestructive) + call util_spill(keeps%rb, discards%rb, lspill_list, ldestructive) + call util_spill(keeps%vb, discards%vb, lspill_list, ldestructive) + call util_spill(keeps%ah, discards%ah, lspill_list, ldestructive) + call util_spill(keeps%aobl, discards%aobl, lspill_list, ldestructive) + call util_spill(keeps%agr, discards%agr, lspill_list, ldestructive) + call util_spill(keeps%atide, discards%atide, lspill_list, ldestructive) + call util_spill(keeps%ir3h, discards%ir3h, lspill_list, ldestructive) + call util_spill(keeps%isperi, discards%isperi, lspill_list, ldestructive) + call util_spill(keeps%peri, discards%peri, lspill_list, ldestructive) + call util_spill(keeps%atp, discards%atp, lspill_list, ldestructive) + call util_spill(keeps%a, discards%a, lspill_list, ldestructive) + call util_spill(keeps%e, discards%e, lspill_list, ldestructive) + call util_spill(keeps%inc, discards%inc, lspill_list, ldestructive) + call util_spill(keeps%capom, discards%capom, lspill_list, ldestructive) + call util_spill(keeps%omega, discards%omega, lspill_list, ldestructive) + call util_spill(keeps%capm, discards%capm, lspill_list, ldestructive) nbody_old = keeps%nbody @@ -4587,24 +3133,24 @@ module subroutine swiftest_util_spill_pl(self, discards, lspill_list, ldestructi select type (discards) ! The standard requires us to select the type of both arguments in order to access all the components class is (swiftest_pl) !> Spill components specific to the massive body class - call swiftest_util_spill(keeps%mass, discards%mass, lspill_list, ldestructive) - call swiftest_util_spill(keeps%Gmass, discards%Gmass, lspill_list, ldestructive) - call swiftest_util_spill(keeps%rhill, discards%rhill, lspill_list, ldestructive) - call swiftest_util_spill(keeps%renc, discards%renc, lspill_list, ldestructive) - call swiftest_util_spill(keeps%radius, discards%radius, lspill_list, ldestructive) - call swiftest_util_spill(keeps%density, discards%density, lspill_list, ldestructive) - call swiftest_util_spill(keeps%rbeg, discards%rbeg, lspill_list, ldestructive) - call swiftest_util_spill(keeps%rend, discards%rend, lspill_list, ldestructive) - call swiftest_util_spill(keeps%vbeg, discards%vbeg, lspill_list, ldestructive) - call swiftest_util_spill(keeps%Ip, discards%Ip, lspill_list, ldestructive) - call swiftest_util_spill(keeps%rot, discards%rot, lspill_list, ldestructive) - call swiftest_util_spill(keeps%k2, discards%k2, lspill_list, ldestructive) - call swiftest_util_spill(keeps%Q, discards%Q, lspill_list, ldestructive) - call swiftest_util_spill(keeps%tlag, discards%tlag, lspill_list, ldestructive) - call swiftest_util_spill(keeps%kin, discards%kin, lspill_list, ldestructive) - call swiftest_util_spill(keeps%lmtiny, discards%lmtiny, lspill_list, ldestructive) - call swiftest_util_spill(keeps%nplenc, discards%nplenc, lspill_list, ldestructive) - call swiftest_util_spill(keeps%ntpenc, discards%ntpenc, lspill_list, ldestructive) + call util_spill(keeps%mass, discards%mass, lspill_list, ldestructive) + call util_spill(keeps%Gmass, discards%Gmass, lspill_list, ldestructive) + call util_spill(keeps%rhill, discards%rhill, lspill_list, ldestructive) + call util_spill(keeps%renc, discards%renc, lspill_list, ldestructive) + call util_spill(keeps%radius, discards%radius, lspill_list, ldestructive) + call util_spill(keeps%density, discards%density, lspill_list, ldestructive) + call util_spill(keeps%rbeg, discards%rbeg, lspill_list, ldestructive) + call util_spill(keeps%rend, discards%rend, lspill_list, ldestructive) + call util_spill(keeps%vbeg, discards%vbeg, lspill_list, ldestructive) + call util_spill(keeps%Ip, discards%Ip, lspill_list, ldestructive) + call util_spill(keeps%rot, discards%rot, lspill_list, ldestructive) + call util_spill(keeps%k2, discards%k2, lspill_list, ldestructive) + call util_spill(keeps%Q, discards%Q, lspill_list, ldestructive) + call util_spill(keeps%tlag, discards%tlag, lspill_list, ldestructive) + call util_spill(keeps%kin, discards%kin, lspill_list, ldestructive) + call util_spill(keeps%lmtiny, discards%lmtiny, lspill_list, ldestructive) + call util_spill(keeps%nplenc, discards%nplenc, lspill_list, ldestructive) + call util_spill(keeps%ntpenc, discards%ntpenc, lspill_list, ldestructive) if (ldestructive .and. allocated(keeps%k_plpl)) deallocate(keeps%k_plpl) @@ -4634,7 +3180,7 @@ module subroutine swiftest_util_spill_tp(self, discards, lspill_list, ldestructi select type(discards) class is (swiftest_tp) !> Spill components specific to the test particle class - call swiftest_util_spill(keeps%nplenc, discards%nplenc, lspill_list, ldestructive) + call util_spill(keeps%nplenc, discards%nplenc, lspill_list, ldestructive) call swiftest_util_spill_body(keeps, discards, lspill_list, ldestructive) class default write(*,*) 'Error! spill method called for incompatible return type on swiftest_tp' @@ -4645,70 +3191,6 @@ module subroutine swiftest_util_spill_tp(self, discards, lspill_list, ldestructi end subroutine swiftest_util_spill_tp - module subroutine swiftest_util_unique_DP(input_array, output_array, index_map) - !! author: David A. Minton - !! - !! Takes an input unsorted integer array and returns a new array of sorted, unique values (DP version) - implicit none - ! Arguments - real(DP), dimension(:), intent(in) :: input_array !! Unsorted input array - real(DP), dimension(:), allocatable, intent(out) :: output_array !! Sorted array of unique values - integer(I4B), dimension(:), allocatable, intent(out) :: index_map !! An array of the same size as input_array that such that any for any index i, output_array(index_map(i)) = input_array(i) - ! Internals - real(DP), dimension(:), allocatable :: unique_array - integer(I4B) :: n - real(DP) :: lo, hi - - allocate(unique_array, mold=input_array) - allocate(index_map(size(input_array))) - lo = minval(input_array) - 1 - hi = maxval(input_array) - - n = 0 - do - n = n + 1 - lo = minval(input_array(:), mask=input_array(:) > lo) - unique_array(n) = lo - where(input_array(:) == lo) index_map(:) = n - if (lo >= hi) exit - enddo - allocate(output_array(n), source=unique_array(1:n)) - - return - end subroutine swiftest_util_unique_DP - - - module subroutine swiftest_util_unique_I4B(input_array, output_array, index_map) - !! author: David A. Minton - !! - !! Takes an input unsorted integer array and returns a new array of sorted, unique values (I4B version) - implicit none - ! Arguments - integer(I4B), dimension(:), intent(in) :: input_array !! Unsorted input array - integer(I4B), dimension(:), allocatable, intent(out) :: output_array !! Sorted array of unique values - integer(I4B), dimension(:), allocatable, intent(out) :: index_map !! An array of the same size as input_array that such that any for any index i, output_array(index_map(i)) = input_array(i) - ! Internals - integer(I4B), dimension(:), allocatable :: unique_array - integer(I4B) :: n, lo, hi - - allocate(unique_array, mold=input_array) - allocate(index_map, mold=input_array) - lo = minval(input_array) - 1 - hi = maxval(input_array) - - n = 0 - do - n = n + 1 - lo = minval(input_array(:), mask=input_array(:) > lo) - unique_array(n) = lo - where(input_array(:) == lo) index_map(:) = n - if (lo >= hi) exit - enddo - allocate(output_array(n), source=unique_array(1:n)) - - return - end subroutine swiftest_util_unique_I4B - module subroutine swiftest_util_valid_id_system(self, param) !! author: David A. Minton @@ -4740,11 +3222,11 @@ module subroutine swiftest_util_valid_id_system(self, param) maxid = maxval(idarr) ! Check to see if the ids are unique - call swiftest_util_unique(idarr, unique_idarr, idmap) + call util_unique(idarr, unique_idarr, idmap) if (size(unique_idarr) == nid) return ! All id values are unique ! Fix any duplicate id values and update the maxid - call swiftest_util_sort(idmap) + call util_sort(idmap) do i = 2, size(idmap) if (idmap(i) == idmap(i-1)) then maxid = maxid + 1 diff --git a/src/symba/symba_util.f90 b/src/symba/symba_util.f90 index ad1d50383..94b6aa9f5 100644 --- a/src/symba/symba_util.f90 +++ b/src/symba/symba_util.f90 @@ -25,12 +25,10 @@ module subroutine symba_util_append_pl(self, source, lsource_mask) select type(source) class is (symba_pl) - associate(nold => self%nbody, nsrc => source%nbody) - call swiftest_util_append(self%levelg, source%levelg, nold, nsrc, lsource_mask) - call swiftest_util_append(self%levelm, source%levelm, nold, nsrc, lsource_mask) + call util_append(self%levelg, source%levelg, lsource_mask=lsource_mask) + call util_append(self%levelm, source%levelm, lsource_mask=lsource_mask) - call swiftest_util_append_pl(self, source, lsource_mask) ! Note: helio_pl does not have its own append method, so we skip back to the base class - end associate + call swiftest_util_append_pl(self, source, lsource_mask) ! Note: helio_pl does not have its own append method, so we skip back to the base class class default write(*,*) "Invalid object passed to the append method. Source must be of class symba_pl or its descendents!" call base_util_exit(FAILURE) @@ -53,12 +51,10 @@ module subroutine symba_util_append_tp(self, source, lsource_mask) select type(source) class is (symba_tp) - associate(nold => self%nbody, nsrc => source%nbody) - call swiftest_util_append(self%levelg, source%levelg, nold, nsrc, lsource_mask) - call swiftest_util_append(self%levelm, source%levelm, nold, nsrc, lsource_mask) + call util_append(self%levelg, source%levelg, lsource_mask=lsource_mask) + call util_append(self%levelm, source%levelm, lsource_mask=lsource_mask) - call swiftest_util_append_tp(self, source, lsource_mask) ! Note: helio_tp does not have its own append method, so we skip back to the base class - end associate + call swiftest_util_append_tp(self, source, lsource_mask) ! Note: helio_tp does not have its own append method, so we skip back to the base class class default write(*,*) "Invalid object passed to the append method. Source must be of class symba_tp or its descendents!" call base_util_exit(FAILURE) @@ -134,8 +130,8 @@ module subroutine symba_util_fill_pl(self, inserts, lfill_list) associate(keeps => self) select type(inserts) class is (symba_pl) - call swiftest_util_fill(keeps%levelg, inserts%levelg, lfill_list) - call swiftest_util_fill(keeps%levelm, inserts%levelm, lfill_list) + call util_fill(keeps%levelg, inserts%levelg, lfill_list) + call util_fill(keeps%levelm, inserts%levelm, lfill_list) call swiftest_util_fill_pl(keeps, inserts, lfill_list) ! Note: helio_pl does not have its own fill method, so we skip back to the base class class default @@ -163,9 +159,9 @@ module subroutine symba_util_fill_tp(self, inserts, lfill_list) associate(keeps => self) select type(inserts) class is (symba_tp) - call swiftest_util_fill(keeps%nplenc, inserts%nplenc, lfill_list) - call swiftest_util_fill(keeps%levelg, inserts%levelg, lfill_list) - call swiftest_util_fill(keeps%levelm, inserts%levelm, lfill_list) + call util_fill(keeps%nplenc, inserts%nplenc, lfill_list) + call util_fill(keeps%levelg, inserts%levelg, lfill_list) + call util_fill(keeps%levelm, inserts%levelm, lfill_list) call swiftest_util_fill_tp(keeps, inserts, lfill_list) ! Note: helio_tp does not have its own fill method, so we skip back to the base class class default @@ -223,8 +219,8 @@ module subroutine symba_util_resize_pl(self, nnew) class(symba_pl), intent(inout) :: self !! SyMBA massive body object integer(I4B), intent(in) :: nnew !! New size neded - call swiftest_util_resize(self%levelg, nnew) - call swiftest_util_resize(self%levelm, nnew) + call util_resize(self%levelg, nnew) + call util_resize(self%levelm, nnew) call swiftest_util_resize_pl(self, nnew) @@ -240,8 +236,8 @@ module subroutine symba_util_resize_tp(self, nnew) class(symba_tp), intent(inout) :: self !! SyMBA test particle object integer(I4B), intent(in) :: nnew !! New size neded - call swiftest_util_resize(self%levelg, nnew) - call swiftest_util_resize(self%levelm, nnew) + call util_resize(self%levelg, nnew) + call util_resize(self%levelm, nnew) call swiftest_util_resize_tp(self, nnew) @@ -423,9 +419,9 @@ module subroutine symba_util_sort_pl(self, sortby, ascending) associate(pl => self, npl => self%nbody) select case(sortby) case("levelg") - call swiftest_util_sort(direction * pl%levelg(1:npl), ind) + call util_sort(direction * pl%levelg(1:npl), ind) case("levelm") - call swiftest_util_sort(direction * pl%levelm(1:npl), ind) + call util_sort(direction * pl%levelm(1:npl), ind) case default ! Look for components in the parent class call swiftest_util_sort_pl(pl, sortby, ascending) @@ -464,11 +460,11 @@ module subroutine symba_util_sort_tp(self, sortby, ascending) associate(tp => self, ntp => self%nbody) select case(sortby) case("nplenc") - call swiftest_util_sort(direction * tp%nplenc(1:ntp), ind) + call util_sort(direction * tp%nplenc(1:ntp), ind) case("levelg") - call swiftest_util_sort(direction * tp%levelg(1:ntp), ind) + call util_sort(direction * tp%levelg(1:ntp), ind) case("levelm") - call swiftest_util_sort(direction * tp%levelm(1:ntp), ind) + call util_sort(direction * tp%levelm(1:ntp), ind) case default ! Look for components in the parent class call swiftest_util_sort_tp(tp, sortby, ascending) return @@ -492,8 +488,8 @@ module subroutine symba_util_sort_rearrange_pl(self, ind) integer(I4B), dimension(:), intent(in) :: ind !! Index array used to restructure the body (should contain all 1:n index values in the desired order) associate(pl => self, npl => self%nbody) - call swiftest_util_sort_rearrange(pl%levelg, ind, npl) - call swiftest_util_sort_rearrange(pl%levelm, ind, npl) + call util_sort_rearrange(pl%levelg, ind, npl) + call util_sort_rearrange(pl%levelm, ind, npl) call swiftest_util_sort_rearrange_pl(pl,ind) end associate @@ -512,9 +508,9 @@ module subroutine symba_util_sort_rearrange_tp(self, ind) integer(I4B), dimension(:), intent(in) :: ind !! Index array used to restructure the body (should contain all 1:n index values in the desired order) associate(tp => self, ntp => self%nbody) - call swiftest_util_sort_rearrange(tp%nplenc, ind, ntp) - call swiftest_util_sort_rearrange(tp%levelg, ind, ntp) - call swiftest_util_sort_rearrange(tp%levelm, ind, ntp) + call util_sort_rearrange(tp%nplenc, ind, ntp) + call util_sort_rearrange(tp%levelg, ind, ntp) + call util_sort_rearrange(tp%levelm, ind, ntp) call swiftest_util_sort_rearrange_tp(tp,ind) end associate @@ -540,8 +536,8 @@ module subroutine symba_util_spill_pl(self, discards, lspill_list, ldestructive) associate(keeps => self) select type(discards) class is (symba_pl) - call swiftest_util_spill(keeps%levelg, discards%levelg, lspill_list, ldestructive) - call swiftest_util_spill(keeps%levelm, discards%levelm, lspill_list, ldestructive) + call util_spill(keeps%levelg, discards%levelg, lspill_list, ldestructive) + call util_spill(keeps%levelm, discards%levelm, lspill_list, ldestructive) call swiftest_util_spill_pl(keeps, discards, lspill_list, ldestructive) class default @@ -571,9 +567,9 @@ module subroutine symba_util_spill_tp(self, discards, lspill_list, ldestructive) associate(keeps => self) select type(discards) class is (symba_tp) - call swiftest_util_spill(keeps%nplenc, discards%nplenc, lspill_list, ldestructive) - call swiftest_util_spill(keeps%levelg, discards%levelg, lspill_list, ldestructive) - call swiftest_util_spill(keeps%levelm, discards%levelm, lspill_list, ldestructive) + call util_spill(keeps%nplenc, discards%nplenc, lspill_list, ldestructive) + call util_spill(keeps%levelg, discards%levelg, lspill_list, ldestructive) + call util_spill(keeps%levelm, discards%levelm, lspill_list, ldestructive) call swiftest_util_spill_tp(keeps, discards, lspill_list, ldestructive) class default diff --git a/src/whm/whm_util.f90 b/src/whm/whm_util.f90 index 13105047e..09fa59d76 100644 --- a/src/whm/whm_util.f90 +++ b/src/whm/whm_util.f90 @@ -24,15 +24,13 @@ module subroutine whm_util_append_pl(self, source, lsource_mask) select type(source) class is (whm_pl) - associate(nold => self%nbody, nsrc => source%nbody) - call swiftest_util_append(self%eta, source%eta, nold, nsrc, lsource_mask) - call swiftest_util_append(self%muj, source%muj, nold, nsrc, lsource_mask) - call swiftest_util_append(self%ir3j, source%ir3j, nold, nsrc, lsource_mask) - call swiftest_util_append(self%xj, source%xj, nold, nsrc, lsource_mask) - call swiftest_util_append(self%vj, source%vj, nold, nsrc, lsource_mask) - - call swiftest_util_append_pl(self, source, lsource_mask) - end associate + call util_append(self%eta, source%eta, lsource_mask=lsource_mask) + call util_append(self%muj, source%muj, lsource_mask=lsource_mask) + call util_append(self%ir3j, source%ir3j, lsource_mask=lsource_mask) + call util_append(self%xj, source%xj, lsource_mask=lsource_mask) + call util_append(self%vj, source%vj, lsource_mask=lsource_mask) + + call swiftest_util_append_pl(self, source, lsource_mask) class default write(*,*) "Invalid object passed to the append method. Source must be of class whm_pl or its descendents" call base_util_exit(FAILURE) @@ -78,11 +76,11 @@ module subroutine whm_util_fill_pl(self, inserts, lfill_list) associate(keeps => self) select type(inserts) class is (whm_pl) - call swiftest_util_fill(keeps%eta, inserts%eta, lfill_list) - call swiftest_util_fill(keeps%muj, inserts%muj, lfill_list) - call swiftest_util_fill(keeps%ir3j, inserts%ir3j, lfill_list) - call swiftest_util_fill(keeps%xj, inserts%xj, lfill_list) - call swiftest_util_fill(keeps%vj, inserts%vj, lfill_list) + call util_fill(keeps%eta, inserts%eta, lfill_list) + call util_fill(keeps%muj, inserts%muj, lfill_list) + call util_fill(keeps%ir3j, inserts%ir3j, lfill_list) + call util_fill(keeps%xj, inserts%xj, lfill_list) + call util_fill(keeps%vj, inserts%vj, lfill_list) call swiftest_util_fill_pl(keeps, inserts, lfill_list) class default @@ -104,11 +102,11 @@ module subroutine whm_util_resize_pl(self, nnew) class(whm_pl), intent(inout) :: self !! WHM massive body object integer(I4B), intent(in) :: nnew !! New size neded - call swiftest_util_resize(self%eta, nnew) - call swiftest_util_resize(self%xj, nnew) - call swiftest_util_resize(self%vj, nnew) - call swiftest_util_resize(self%muj, nnew) - call swiftest_util_resize(self%ir3j, nnew) + call util_resize(self%eta, nnew) + call util_resize(self%xj, nnew) + call util_resize(self%vj, nnew) + call util_resize(self%muj, nnew) + call util_resize(self%ir3j, nnew) call swiftest_util_resize_pl(self, nnew) @@ -255,11 +253,11 @@ module subroutine whm_util_sort_pl(self, sortby, ascending) associate(pl => self, npl => self%nbody) select case(sortby) case("eta") - call swiftest_util_sort(direction * pl%eta(1:npl), ind) + call util_sort(direction * pl%eta(1:npl), ind) case("muj") - call swiftest_util_sort(direction * pl%muj(1:npl), ind) + call util_sort(direction * pl%muj(1:npl), ind) case("ir3j") - call swiftest_util_sort(direction * pl%ir3j(1:npl), ind) + call util_sort(direction * pl%ir3j(1:npl), ind) case("xj", "vj") write(*,*) 'Cannot sort by ' // trim(adjustl(sortby)) // '. Component not sortable!' case default @@ -287,11 +285,11 @@ module subroutine whm_util_sort_rearrange_pl(self, ind) if (self%nbody == 0) return associate(pl => self, npl => self%nbody) - call swiftest_util_sort_rearrange(pl%eta, ind, npl) - call swiftest_util_sort_rearrange(pl%xj, ind, npl) - call swiftest_util_sort_rearrange(pl%vj, ind, npl) - call swiftest_util_sort_rearrange(pl%muj, ind, npl) - call swiftest_util_sort_rearrange(pl%ir3j, ind, npl) + call util_sort_rearrange(pl%eta, ind, npl) + call util_sort_rearrange(pl%xj, ind, npl) + call util_sort_rearrange(pl%vj, ind, npl) + call util_sort_rearrange(pl%muj, ind, npl) + call util_sort_rearrange(pl%ir3j, ind, npl) call swiftest_util_sort_rearrange_pl(pl,ind) end associate @@ -316,11 +314,11 @@ module subroutine whm_util_spill_pl(self, discards, lspill_list, ldestructive) associate(keeps => self) select type(discards) class is (whm_pl) - call swiftest_util_spill(keeps%eta, discards%eta, lspill_list, ldestructive) - call swiftest_util_spill(keeps%muj, discards%muj, lspill_list, ldestructive) - call swiftest_util_spill(keeps%ir3j, discards%ir3j, lspill_list, ldestructive) - call swiftest_util_spill(keeps%xj, discards%xj, lspill_list, ldestructive) - call swiftest_util_spill(keeps%vj, discards%vj, lspill_list, ldestructive) + call util_spill(keeps%eta, discards%eta, lspill_list, ldestructive) + call util_spill(keeps%muj, discards%muj, lspill_list, ldestructive) + call util_spill(keeps%ir3j, discards%ir3j, lspill_list, ldestructive) + call util_spill(keeps%xj, discards%xj, lspill_list, ldestructive) + call util_spill(keeps%vj, discards%vj, lspill_list, ldestructive) call swiftest_util_spill_pl(keeps, discards, lspill_list, ldestructive) class default