diff --git a/src/modules/swiftest_classes.f90 b/src/modules/swiftest_classes.f90 index 3711e5295..c8d05850f 100644 --- a/src/modules/swiftest_classes.f90 +++ b/src/modules/swiftest_classes.f90 @@ -85,7 +85,7 @@ module swiftest_classes !******************************************************************************************************************************** !> A concrete lass for the central body in a Swiftest simulation type, abstract, extends(swiftest_base) :: swiftest_cb - character(len=STRMAX) :: name !! Non-unique name + character(len=STRMAX) :: name !! Non-unique name integer(I4B) :: id = 0 !! External identifier (unique) real(DP) :: mass = 0.0_DP !! Central body mass (units MU) real(DP) :: Gmass = 0.0_DP !! Central mass gravitational term G * mass (units GU * MU) @@ -801,7 +801,46 @@ module subroutine util_copy_fill_tp(self, inserts, lfill_list) class(swiftest_body), intent(in) :: inserts !! Swiftest body object to be inserted logical, dimension(:), intent(in) :: lfill_list !! Logical array of bodies to merge into the keeps end subroutine util_copy_fill_tp + end interface + + interface util_fill + module subroutine 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 arrays to insert into keep + logical, dimension(:), intent(in) :: lfill_list !! Logical array of bodies to merge into the keeps + end subroutine util_fill_arr_char_string + + module subroutine 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 arrays to insert into keep + logical, dimension(:), intent(in) :: lfill_list !! Logical array of bodies to merge into the keeps + end subroutine util_fill_arr_DP + + module subroutine 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 arrays to insert into keep + logical, dimension(:), intent(in) :: lfill_list !! Logical array of bodies to merge into the keeps + end subroutine util_fill_arr_DPvec + module subroutine 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 arrays to insert into keep + logical, dimension(:), intent(in) :: lfill_list !! Logical array of bodies to merge into the keeps + end subroutine util_fill_arr_I4B + + module subroutine 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 arrays to insert into keep + logical, dimension(:), intent(in) :: lfill_list !! Logical array of bodies to merge into the keeps + end subroutine util_fill_arr_logical + end interface + + interface module subroutine util_peri_tp(self, system, param) implicit none class(swiftest_tp), intent(inout) :: self !! Swiftest test particle object diff --git a/src/rmvs/rmvs_util.f90 b/src/rmvs/rmvs_util.f90 index 863033c17..90949e593 100644 --- a/src/rmvs/rmvs_util.f90 +++ b/src/rmvs/rmvs_util.f90 @@ -20,9 +20,10 @@ module subroutine rmvs_util_copy_fill_pl(self, inserts, lfill_list) select type(inserts) class is (rmvs_pl) - keeps%nenc(:) = unpack(keeps%nenc(:), .not.lfill_list(:), keeps%nenc(:)) - keeps%nenc(:) = unpack(inserts%nenc(:), lfill_list(:), keeps%nenc(:)) - + 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) + call whm_util_copy_fill_pl(keeps, inserts, lfill_list) class default write(*,*) 'Error! spill method called for incompatible return type on rmvs_pl' @@ -49,14 +50,9 @@ module subroutine rmvs_util_copy_fill_tp(self, inserts, lfill_list) select type(inserts) class is (rmvs_tp) - keeps%lperi(:) = unpack(keeps%lperi(:), .not.lfill_list(:), keeps%lperi(:)) - keeps%lperi(:) = unpack(inserts%lperi(:), lfill_list(:), keeps%lperi(:)) - - keeps%plperP(:) = unpack(keeps%plperP(:), .not.lfill_list(:), keeps%plperP(:)) - keeps%plperP(:) = unpack(inserts%plperP(:), lfill_list(:), keeps%plperP(:)) - - keeps%plencP(:) = unpack(keeps%plencP(:), .not.lfill_list(:), keeps%plencP(:)) - keeps%plencP(:) = unpack(inserts%plencP(:), lfill_list(:), keeps%plencP(:)) + 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 util_copy_fill_tp(keeps, inserts, lfill_list) class default diff --git a/src/symba/symba_util.f90 b/src/symba/symba_util.f90 index 7d65665f8..96781555c 100644 --- a/src/symba/symba_util.f90 +++ b/src/symba/symba_util.f90 @@ -17,35 +17,16 @@ module subroutine symba_util_copy_fill_pl(self, inserts, lfill_list) associate(keeps => self) select type(inserts) class is (symba_pl) - keeps%lcollision(:) = unpack(keeps%lcollision(:), .not.lfill_list(:), keeps%lcollision(:)) - keeps%lcollision(:) = unpack(inserts%lcollision(:), lfill_list(:), keeps%lcollision(:)) - - keeps%lencounter(:) = unpack(keeps%lencounter(:), .not.lfill_list(:), keeps%lencounter(:)) - keeps%lencounter(:) = unpack(inserts%lencounter(:), lfill_list(:), keeps%lencounter(:)) - - keeps%lmtiny(:) = unpack(keeps%lmtiny(:), .not.lfill_list(:), keeps%lmtiny(:)) - keeps%lmtiny(:) = unpack(inserts%lmtiny(:), lfill_list(:), keeps%lmtiny(:)) - - keeps%nplenc(:) = unpack(keeps%nplenc(:), .not.lfill_list(:), keeps%nplenc(:)) - keeps%nplenc(:) = unpack(inserts%nplenc(:), lfill_list(:), keeps%nplenc(:)) - - keeps%nplenc(:) = unpack(keeps%nplenc(:), .not.lfill_list(:), keeps%nplenc(:)) - keeps%ntpenc(:) = unpack(inserts%ntpenc(:), lfill_list(:), keeps%ntpenc(:)) - - keeps%levelg(:) = unpack(keeps%levelg(:), .not.lfill_list(:), keeps%levelg(:)) - keeps%levelg(:) = unpack(inserts%levelg(:), lfill_list(:), keeps%levelg(:)) - - keeps%levelm(:) = unpack(keeps%levelm(:), .not.lfill_list(:), keeps%levelm(:)) - keeps%levelm(:) = unpack(inserts%levelm(:), lfill_list(:), keeps%levelm(:)) - - keeps%isperi(:) = unpack(keeps%isperi(:), .not.lfill_list(:), keeps%isperi(:)) - keeps%isperi(:) = unpack(inserts%isperi(:), lfill_list(:), keeps%isperi(:)) - - keeps%peri(:) = unpack(keeps%peri(:), .not.lfill_list(:), keeps%peri(:)) - keeps%peri(:) = unpack(inserts%peri(:), lfill_list(:), keeps%peri(:)) - - keeps%atp(:) = unpack(keeps%atp(:), .not.lfill_list(:), keeps%atp(:)) - keeps%atp(:) = unpack(inserts%atp(:), lfill_list(:), keeps%atp(:)) + call util_fill(keeps%lcollision, inserts%lcollision, lfill_list) + call util_fill(keeps%lencounter, inserts%lencounter, lfill_list) + call util_fill(keeps%lmtiny, inserts%lmtiny, lfill_list) + call util_fill(keeps%nplenc, inserts%nplenc, lfill_list) + call util_fill(keeps%ntpenc, inserts%ntpenc, lfill_list) + call util_fill(keeps%levelg, inserts%levelg, lfill_list) + call util_fill(keeps%levelm, inserts%levelm, 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) keeps%kin(:) = unpack(keeps%kin(:), .not.lfill_list(:), keeps%kin(:)) keeps%kin(:) = unpack(inserts%kin(:), lfill_list(:), keeps%kin(:)) @@ -77,14 +58,9 @@ module subroutine symba_util_copy_fill_tp(self, inserts, lfill_list) associate(keeps => self) select type(inserts) class is (symba_tp) - keeps%nplenc(:) = unpack(keeps%nplenc(:), .not.lfill_list(:), keeps%nplenc(:)) - keeps%nplenc(:) = unpack(inserts%nplenc(:), lfill_list(:), keeps%nplenc(:)) - - keeps%levelg(:) = unpack(keeps%levelg(:), .not.lfill_list(:), keeps%levelg(:)) - keeps%levelg(:) = unpack(inserts%levelg(:), lfill_list(:), keeps%levelg(:)) - - keeps%levelm(:) = unpack(keeps%levelm(:), .not.lfill_list(:), keeps%levelm(:)) - keeps%levelm(:) = unpack(inserts%levelm(:), lfill_list(:), keeps%levelm(:)) + 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 util_copy_fill_tp(keeps, inserts, lfill_list) class default diff --git a/src/util/util_copy.f90 b/src/util/util_copy.f90 index 261a490fd..60cf568b4 100644 --- a/src/util/util_copy.f90 +++ b/src/util/util_copy.f90 @@ -2,249 +2,6 @@ use swiftest contains - module subroutine util_copy_fill_body(self, inserts, lfill_list) - !! author: David A. Minton - !! - !! Insert new Swiftest generic particle structure into an old one. - !! This is the inverse of a fill operation. - implicit none - ! Arguments - class(swiftest_body), intent(inout) :: self !! Swiftest generic body object - class(swiftest_body), intent(in) :: inserts !! Inserted object - logical, dimension(:), intent(in) :: lfill_list !! Logical array of bodies to merge into the keeps - ! internals - integer(I4B) :: i - - ! For each component, pack the discarded bodies into the discard object and do the inverse with the keeps - !> Spill all the common components - associate(keeps => self) - keeps%id(:) = unpack(keeps%id(:), .not.lfill_list(:), keeps%id(:)) - keeps%id(:) = unpack(inserts%id(:), lfill_list(:), keeps%id(:)) - - keeps%name(:) = unpack(keeps%name(:), .not.lfill_list(:), keeps%name(:)) - keeps%name(:) = unpack(inserts%name(:), lfill_list(:), keeps%name(:)) - - keeps%status(:) = unpack(keeps%status(:), .not.lfill_list(:), keeps%status(:)) - keeps%status(:) = unpack(inserts%status(:), lfill_list(:), keeps%status(:)) - - keeps%ldiscard(:) = unpack(keeps%ldiscard(:), .not.lfill_list(:), keeps%ldiscard(:)) - keeps%ldiscard(:) = unpack(inserts%ldiscard(:), lfill_list(:), keeps%ldiscard(:)) - - keeps%mu(:) = unpack(keeps%mu(:), .not.lfill_list(:), keeps%mu(:)) - keeps%mu(:) = unpack(inserts%mu(:), lfill_list(:), keeps%mu(:)) - - keeps%lmask(:) = unpack(keeps%lmask(:), .not.lfill_list(:), keeps%ldiscard(:)) - keeps%lmask(:) = unpack(inserts%lmask(:), lfill_list(:), keeps%ldiscard(:)) - - do i = 1, NDIM - keeps%xh(i, :) = unpack(keeps%xh(i, :), .not.lfill_list(:), keeps%xh(i, :)) - keeps%xh(i, :) = unpack(inserts%xh(i, :), lfill_list(:), keeps%xh(i, :)) - - keeps%vh(i, :) = unpack(keeps%vh(i, :), .not.lfill_list(:), keeps%vh(i, :)) - keeps%vh(i, :) = unpack(inserts%vh(i, :), lfill_list(:), keeps%vh(i, :)) - - keeps%xb(i, :) = unpack(keeps%xb(i, :), .not.lfill_list(:), keeps%xb(i, :)) - keeps%xb(i, :) = unpack(inserts%xb(i, :), lfill_list(:), keeps%xb(i, :)) - - keeps%vb(i, :) = unpack(keeps%vb(i, :), .not.lfill_list(:), keeps%vb(i, :)) - keeps%vb(i, :) = unpack(inserts%vb(i, :), lfill_list(:), keeps%vb(i, :)) - - keeps%ah(i, :) = unpack(keeps%ah(i, :), .not.lfill_list(:), keeps%ah(i, :)) - keeps%ah(i, :) = unpack(inserts%ah(i, :), lfill_list(:), keeps%ah(i, :)) - end do - - if (allocated(keeps%aobl)) then - do i = 1, NDIM - keeps%aobl(i, :) = unpack(keeps%aobl(i, :), .not.lfill_list(:), keeps%aobl(i, :)) - keeps%aobl(i, :) = unpack(inserts%aobl(i, :), lfill_list(:), keeps%aobl(i, :)) - end do - end if - - if (allocated(keeps%agr)) then - do i = 1, NDIM - keeps%agr(i, :) = unpack(keeps%agr(i, :), .not.lfill_list(:), keeps%agr(i, :)) - keeps%agr(i, :) = unpack(inserts%agr(i, :), lfill_list(:), keeps%agr(i, :)) - end do - end if - - if (allocated(keeps%atide)) then - do i = 1, NDIM - keeps%atide(i, :) = unpack(keeps%atide(i, :), .not.lfill_list(:), keeps%atide(i, :)) - keeps%atide(i, :) = unpack(inserts%atide(i, :), lfill_list(:), keeps%atide(i, :)) - end do - end if - - if (allocated(keeps%a)) then - keeps%a(:) = unpack(keeps%a(:), .not.lfill_list(:), keeps%a(:)) - keeps%a(:) = unpack(inserts%a(:), lfill_list(:), keeps%a(:)) - end if - - if (allocated(keeps%e)) then - keeps%e(:) = unpack(keeps%e(:), .not.lfill_list(:), keeps%e(:)) - keeps%e(:) = unpack(inserts%e(:), lfill_list(:), keeps%e(:)) - end if - - if (allocated(keeps%inc)) then - keeps%inc(:) = unpack(keeps%inc(:), .not.lfill_list(:), keeps%inc(:)) - keeps%inc(:) = unpack(inserts%inc(:), lfill_list(:), keeps%inc(:)) - end if - - if (allocated(keeps%capom)) then - keeps%capom(:) = unpack(keeps%capom(:),.not.lfill_list(:), keeps%capom(:)) - keeps%capom(:) = unpack(inserts%capom(:),lfill_list(:), keeps%capom(:)) - end if - - if (allocated(keeps%omega)) then - keeps%omega(:) = unpack(keeps%omega(:),.not.lfill_list(:), keeps%omega(:)) - keeps%omega(:) = unpack(inserts%omega(:),lfill_list(:), keeps%omega(:)) - end if - - if (allocated(keeps%capm)) then - keeps%capm(:) = unpack(keeps%capm(:), .not.lfill_list(:), keeps%capm(:)) - keeps%capm(:) = unpack(inserts%capm(:), lfill_list(:), keeps%capm(:)) - end if - - ! This is the base class, so will be the last to be called in the cascade. - keeps%nbody = size(keeps%id(:)) - end associate - - return - end subroutine util_copy_fill_body - - - module subroutine util_copy_fill_pl(self, inserts, lfill_list) - !! author: David A. Minton - !! - !! Insert new Swiftest massive body structure into an old one. - !! This is the inverse of a fill operation. - implicit none - ! Arguments - class(swiftest_pl), intent(inout) :: self !! Swiftest massive body object - class(swiftest_body), intent(in) :: inserts !! Swiftest body object to be inserted - logical, dimension(:), intent(in) :: lfill_list !! Logical array of bodies to merge into the keeps - ! Internals - integer(I4B) :: i - - associate(keeps => self) - - 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) - !> Spill components specific to the massive body class - keeps%mass(:) = unpack(keeps%mass(:),.not.lfill_list(:), keeps%mass(:)) - keeps%mass(:) = unpack(inserts%mass(:),lfill_list(:), keeps%mass(:)) - - keeps%Gmass(:) = unpack(keeps%Gmass(:),.not.lfill_list(:), keeps%Gmass(:)) - keeps%Gmass(:) = unpack(inserts%Gmass(:),lfill_list(:), keeps%Gmass(:)) - - keeps%rhill(:) = unpack(keeps%rhill(:),.not.lfill_list(:), keeps%rhill(:)) - keeps%rhill(:) = unpack(inserts%rhill(:),lfill_list(:), keeps%rhill(:)) - - if (allocated(keeps%radius) .and. allocated(inserts%radius)) then - keeps%radius(:) = unpack(keeps%radius(:),.not.lfill_list(:), keeps%radius(:)) - keeps%radius(:) = unpack(inserts%radius(:),lfill_list(:), keeps%radius(:)) - end if - - if (allocated(keeps%density) .and. allocated(inserts%density)) then - keeps%density(:) = unpack(keeps%density(:),.not.lfill_list(:), keeps%density(:)) - keeps%density(:) = unpack(inserts%density(:),lfill_list(:), keeps%density(:)) - end if - - if (allocated(keeps%k2) .and. allocated(inserts%k2)) then - keeps%k2(:) = unpack(keeps%k2(:),.not.lfill_list(:), keeps%k2(:)) - keeps%k2(:) = unpack(inserts%k2(:),lfill_list(:), keeps%k2(:)) - end if - - if (allocated(keeps%Q) .and. allocated(inserts%Q)) then - keeps%Q(:) = unpack(keeps%Q(:),.not.lfill_list(:), keeps%Q(:)) - keeps%Q(:) = unpack(inserts%Q(:),lfill_list(:), keeps%Q(:)) - end if - - if (allocated(keeps%tlag) .and. allocated(inserts%tlag)) then - keeps%tlag(:) = unpack(keeps%tlag(:),.not.lfill_list(:), keeps%tlag(:)) - keeps%tlag(:) = unpack(inserts%tlag(:),lfill_list(:), keeps%tlag(:)) - end if - - if (allocated(keeps%xbeg) .and. allocated(inserts%xbeg)) then - do i = 1, NDIM - keeps%xbeg(i, :) = unpack(keeps%xbeg(i, :), .not.lfill_list(:), keeps%xbeg(i, :)) - keeps%xbeg(i, :) = unpack(inserts%xbeg(i, :), lfill_list(:), keeps%xbeg(i, :)) - end do - end if - - if (allocated(keeps%xend) .and. allocated(inserts%xend)) then - do i = 1, NDIM - keeps%xend(i, :) = unpack(keeps%xend(i, :), .not.lfill_list(:), keeps%xend(i, :)) - keeps%xend(i, :) = unpack(inserts%xend(i, :), lfill_list(:), keeps%xend(i, :)) - end do - end if - - if (allocated(keeps%vbeg) .and. allocated(inserts%vbeg)) then - do i = 1, NDIM - keeps%vbeg(i, :) = unpack(keeps%vbeg(i, :), .not.lfill_list(:), keeps%vbeg(i, :)) - keeps%vbeg(i, :) = unpack(inserts%vbeg(i, :), lfill_list(:), keeps%vbeg(i, :)) - end do - end if - - if (allocated(keeps%Ip) .and. allocated(inserts%Ip)) then - do i = 1, NDIM - keeps%Ip(i, :) = unpack(keeps%Ip(i, :), .not.lfill_list(:), keeps%Ip(i, :)) - keeps%Ip(i, :) = unpack(inserts%Ip(i, :), lfill_list(:), keeps%Ip(i, :)) - end do - end if - - if (allocated(keeps%rot) .and. allocated(inserts%rot)) then - do i = 1, NDIM - keeps%rot(i, :) = unpack(keeps%rot(i, :), .not.lfill_list(:), keeps%rot(i, :)) - keeps%rot(i, :) = unpack(inserts%rot(i, :), lfill_list(:), keeps%rot(i, :)) - end do - end if - - keeps%ldiscard(:) = unpack(inserts%ldiscard(:), lfill_list(:), keeps%ldiscard(:)) - - call util_copy_fill_body(keeps, inserts, lfill_list) - class default - write(*,*) 'Error! fill method called for incompatible return type on swiftest_pl' - end select - end associate - - return - end subroutine util_copy_fill_pl - - - module subroutine util_copy_fill_tp(self, inserts, lfill_list) - !! author: David A. Minton - !! - !! Insert new Swiftest test particle structure into an old one. - !! This is the inverse of a fill operation. - implicit none - ! Arguments - class(swiftest_tp), intent(inout) :: self !! Swiftest test particle object - class(swiftest_body), intent(in) :: inserts !! Swiftest body object to be inserted - logical, dimension(:), intent(in) :: lfill_list !! Logical array of bodies to merge into the keeps - - associate(keeps => self) - select type(inserts) - class is (swiftest_tp) - !> Spill components specific to the test particle class - keeps%isperi(:) = unpack(keeps%isperi(:), .not.lfill_list(:), keeps%isperi(:)) - keeps%isperi(:) = unpack(inserts%isperi(:), lfill_list(:), keeps%isperi(:)) - - keeps%peri(:) = unpack(keeps%peri(:), .not.lfill_list(:), keeps%peri(:)) - keeps%peri(:) = unpack(inserts%peri(:), lfill_list(:), keeps%peri(:)) - - keeps%atp(:) = unpack(keeps%atp(:), .not.lfill_list(:), keeps%atp(:)) - keeps%atp(:) = unpack(inserts%atp(:), lfill_list(:), keeps%atp(:)) - - call util_copy_fill_body(keeps, inserts, lfill_list) - class default - write(*,*) 'Error! fill method called for incompatible return type on swiftest_tp' - end select - end associate - - return - end subroutine util_copy_fill_tp - - module subroutine util_copy_into_body(self, source, param, lsource_mask) !! author: David A. Minton !! diff --git a/src/util/util_fill.f90 b/src/util/util_fill.f90 new file mode 100644 index 000000000..47981c038 --- /dev/null +++ b/src/util/util_fill.f90 @@ -0,0 +1,218 @@ +submodule (swiftest_classes) s_util_fill + use swiftest +contains + + module subroutine 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 arrays 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 util_fill_arr_char_string + + module subroutine 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 arrays 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 util_fill_arr_DP + + module subroutine 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 arrays 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 util_fill_arr_DPvec + + module subroutine 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 arrays 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 util_fill_arr_I4B + + module subroutine 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 arrays 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 util_fill_arr_logical + + + module subroutine util_copy_fill_body(self, inserts, lfill_list) + !! author: David A. Minton + !! + !! Insert new Swiftest generic particle structure into an old one. + !! This is the inverse of a spill operation. + implicit none + ! Arguments + class(swiftest_body), intent(inout) :: self !! Swiftest generic body object + class(swiftest_body), intent(in) :: inserts !! Inserted object + logical, dimension(:), intent(in) :: lfill_list !! Logical array of bodies to merge into the keeps + ! internals + integer(I4B) :: i + + ! 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 util_fill(keeps%id, inserts%id, lfill_list) + call util_fill(keeps%name, inserts%name, lfill_list) + call util_fill(keeps%status, inserts%status, lfill_list) + call util_fill(keeps%ldiscard, inserts%ldiscard, lfill_list) + call util_fill(keeps%mu, inserts%mu, lfill_list) + call util_fill(keeps%xh, inserts%xh, lfill_list) + call util_fill(keeps%vh, inserts%vh, lfill_list) + call util_fill(keeps%xb, inserts%xb, 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%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(:)) + end associate + + return + end subroutine util_copy_fill_body + + + module subroutine util_copy_fill_pl(self, inserts, lfill_list) + !! author: David A. Minton + !! + !! Insert new Swiftest massive body structure into an old one. + !! This is the inverse of a spill operation. + implicit none + ! Arguments + class(swiftest_pl), intent(inout) :: self !! Swiftest massive body object + class(swiftest_body), intent(in) :: inserts !! Swiftest body object to be inserted + logical, dimension(:), intent(in) :: lfill_list !! Logical array of bodies to merge into the keeps + ! Internals + integer(I4B) :: i + + associate(keeps => self) + + 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 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%radius, inserts%radius, lfill_list) + call util_fill(keeps%density, inserts%density, 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%xbeg, inserts%xbeg, 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_copy_fill_body(keeps, inserts, lfill_list) + class default + write(*,*) 'Error! fill method called for incompatible return type on swiftest_pl' + end select + end associate + + return + end subroutine util_copy_fill_pl + + + module subroutine util_copy_fill_tp(self, inserts, lfill_list) + !! author: David A. Minton + !! + !! Insert new Swiftest test particle structure into an old one. + !! This is the inverse of a fill operation. + implicit none + ! Arguments + class(swiftest_tp), intent(inout) :: self !! Swiftest test particle object + class(swiftest_body), intent(in) :: inserts !! Swiftest body object to be inserted + logical, dimension(:), intent(in) :: lfill_list !! Logical array of bodies to merge into the keeps + + associate(keeps => self) + select type(inserts) + class is (swiftest_tp) + !> Spill components specific to the test particle class + 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_copy_fill_body(keeps, inserts, lfill_list) + class default + write(*,*) 'Error! fill method called for incompatible return type on swiftest_tp' + end select + end associate + + return + end subroutine util_copy_fill_tp + +end submodule s_util_fill \ No newline at end of file diff --git a/src/whm/whm_util.f90 b/src/whm/whm_util.f90 index aaad01e84..ab5461f2d 100644 --- a/src/whm/whm_util.f90 +++ b/src/whm/whm_util.f90 @@ -63,23 +63,12 @@ module subroutine whm_util_copy_fill_pl(self, inserts, lfill_list) associate(keeps => self) select type(inserts) class is (whm_pl) - keeps%eta(:) = unpack(keeps%eta(:), .not.lfill_list(:), keeps%eta(:)) - keeps%eta(:) = unpack(inserts%eta(:), lfill_list(:), keeps%eta(:)) - - keeps%muj(:) = unpack(keeps%muj(:), .not.lfill_list(:), keeps%muj(:)) - keeps%muj(:) = unpack(inserts%muj(:), lfill_list(:), keeps%muj(:)) - - keeps%ir3j(:) = unpack(keeps%ir3j(:), .not.lfill_list(:), keeps%ir3j(:)) - keeps%ir3j(:) = unpack(inserts%ir3j(:), lfill_list(:), keeps%ir3j(:)) - - - do i = 1, NDIM - keeps%xj(i, :) = unpack(keeps%xj(i, :), .not.lfill_list(:), keeps%xj(i, :)) - keeps%xj(i, :) = unpack(inserts%xj(i, :), lfill_list(:), keeps%xj(i, :)) - - keeps%vj(i, :) = unpack(keeps%vj(i, :), .not.lfill_list(:), keeps%vj(i, :)) - keeps%vj(i, :) = unpack(inserts%vj(i, :), lfill_list(:), keeps%vj(i, :)) - end do + 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 util_copy_fill_pl(keeps, inserts, lfill_list) class default write(*,*) 'Error! fill method called for incompatible return type on whm_pl'