diff --git a/src/modules/symba_classes.f90 b/src/modules/symba_classes.f90 index ca51def4c..75d8e0e73 100644 --- a/src/modules/symba_classes.f90 +++ b/src/modules/symba_classes.f90 @@ -78,6 +78,7 @@ module symba_classes procedure :: fill => symba_util_fill_pl !! "Fills" bodies from one object into another depending on the results of a mask (uses the UNPACK intrinsic) procedure :: get_peri => symba_util_peri_pl !! Determine system pericenter passages for massive bodies procedure :: rearray => symba_util_rearray_pl !! Clean up the massive body structures to remove discarded bodies and add new bodies + procedure :: reset_kinship => symba_util_reset_kinship !! Resets the kinship status of bodies procedure :: resize => symba_util_resize_pl !! Checks the current size of a SyMBA massive body against the requested size and resizes it if it is too small. procedure :: sort => symba_util_sort_pl !! Sorts body arrays by a sortable componen procedure :: rearrange => symba_util_sort_rearrange_pl !! Rearranges the order of array elements of body based on an input index array. Used in sorting methods @@ -573,6 +574,14 @@ module subroutine symba_util_rearray_pl(self, system, param) class(symba_nbody_system), intent(inout) :: system !! SyMBA nbody system object class(symba_parameters), intent(inout) :: param !! Current run configuration parameters with SyMBA additions end subroutine symba_util_rearray_pl + + module subroutine symba_util_reset_kinship(self, idx) + implicit none + class(symba_pl), intent(inout) :: self !! SyMBA massive body object + integer(I4B), dimension(:), intent(in) :: idx !! Index array of bodies to reset + integer(I4B) :: i, j + end subroutine symba_util_reset_kinship + end interface interface util_resize diff --git a/src/symba/symba_collision.f90 b/src/symba/symba_collision.f90 index 7a2a38023..ea93948f3 100644 --- a/src/symba/symba_collision.f90 +++ b/src/symba/symba_collision.f90 @@ -90,6 +90,7 @@ module function symba_collision_casedisruption(system, param, family, x, v, mass pl%status(family(:)) = status pl%ldiscard(family(:)) = .false. pl%lcollision(family(:)) = .false. + call pl%reset_kinship(family(:)) end select else ! Populate the list of new bodies @@ -201,12 +202,7 @@ module function symba_collision_casehitandrun(system, param, family, x, v, mass, pl%status(family(:)) = ACTIVE pl%ldiscard(family(:)) = .false. pl%lcollision(family(:)) = .false. - pl%kin(family(:))%parent = family(:) - pl%kin(family(:))%nchild = 0 - do j = 1, size(family(:)) - i = family(j) - if (allocated(pl%kin(i)%child)) deallocate(pl%kin(i)%child) - end do + call pl%reset_kinship(family(:)) end select else status = HIT_AND_RUN_DISRUPT @@ -414,12 +410,7 @@ module function symba_collision_casesupercatastrophic(system, param, family, x, pl%status(family(:)) = status pl%ldiscard(family(:)) = .false. pl%lcollision(family(:)) = .false. - pl%kin(family(:))%parent = family(:) - pl%kin(family(:))%nchild = 0 - do j = 1, size(family(:)) - i = family(j) - if (allocated(pl%kin(i)%child)) deallocate(pl%kin(i)%child) - end do + call pl%reset_kinship(family(:)) end select else ! Populate the list of new bodies @@ -874,8 +865,7 @@ module subroutine symba_collision_make_family_pl(self, idx) pl%kin(j)%parent = index_parent end do end if - if (allocated(pl%kin(index_child)%child)) deallocate(pl%kin(index_child)%child) - pl%kin(index_child)%nchild = 0 + call pl%reset_kinship([index_child]) ! Add the new child to its parent pl%kin(index_child)%parent = index_parent temp(nchild_new) = index_child @@ -906,7 +896,7 @@ subroutine symba_collision_mergeaddsub(system, param, family, id_frag, Ip_frag, ! Internals integer(I4B) :: i, ibiggest, ismallest, iother, nstart, nend, nfamily, nfrag logical, dimension(system%pl%nbody) :: lmask - class(symba_pl), allocatable :: plnew + class(symba_pl), allocatable :: plnew, plsub character(*), parameter :: FRAGFMT = '("Newbody",I0.7)' character(len=NAMELEN) :: newname @@ -1022,15 +1012,18 @@ subroutine symba_collision_mergeaddsub(system, param, family, id_frag, Ip_frag, pl%status(family(:)) = MERGED pl%ldiscard(family(:)) = .true. pl%lcollision(family(:)) = .true. + call pl%reset_kinship(family(:)) lmask(:) = .false. lmask(family(:)) = .true. call plnew%setup(0, param) - call pl%spill(plnew, lmask, ldestructive=.false.) + + allocate(plsub, mold=pl) + call pl%spill(plsub, lmask, ldestructive=.false.) nstart = pl_discards%nbody + 1 nend = pl_discards%nbody + nfamily - call pl_discards%append(plnew, lsource_mask=[(.true., i = 1, nfamily)]) + call pl_discards%append(plsub, lsource_mask=[(.true., i = 1, nfamily)]) ! Record how many bodies were subtracted in this event pl_discards%ncomp(nstart:nend) = nfamily @@ -1042,7 +1035,7 @@ subroutine symba_collision_mergeaddsub(system, param, family, id_frag, Ip_frag, return end subroutine symba_collision_mergeaddsub - + module subroutine symba_collision_resolve_fragmentations(self, system, param) !! author: David A. Minton diff --git a/src/symba/symba_setup.f90 b/src/symba/symba_setup.f90 index 9c1fdb343..3a0543ee1 100644 --- a/src/symba/symba_setup.f90 +++ b/src/symba/symba_setup.f90 @@ -110,8 +110,7 @@ module subroutine symba_setup_pl(self, n, param) self%isperi(:) = 0 self%peri(:) = 0.0_DP self%atp(:) = 0.0_DP - self%kin(:)%nchild = 0 - self%kin(:)%parent = [(i, i=1, n)] + call self%reset_kinship([(i, i=1, n)]) return end subroutine symba_setup_pl diff --git a/src/symba/symba_step.f90 b/src/symba/symba_step.f90 index fe2d74b1a..7ab8113c0 100644 --- a/src/symba/symba_step.f90 +++ b/src/symba/symba_step.f90 @@ -236,11 +236,7 @@ module subroutine symba_step_reset_system(self, param) associate(npl => pl%nbody, ntp => tp%nbody) if (npl > 0) then pl%lcollision(1:npl) = .false. - pl%kin(1:npl)%parent = [(i, i=1, npl)] - pl%kin(1:npl)%nchild = 0 - do i = 1, npl - if (allocated(pl%kin(i)%child)) deallocate(pl%kin(i)%child) - end do + call pl%reset_kinship([(i, i=1, npl)]) pl%nplenc(1:npl) = 0 pl%ntpenc(1:npl) = 0 pl%levelg(1:npl) = -1 diff --git a/src/symba/symba_util.f90 b/src/symba/symba_util.f90 index ccb9b37ec..2ff7dbf18 100644 --- a/src/symba/symba_util.f90 +++ b/src/symba/symba_util.f90 @@ -476,11 +476,7 @@ module subroutine symba_util_rearray_pl(self, system, param) call pl%index(param) ! Reset the kinship trackers - pl%kin(1:npl)%nchild = 0 - pl%kin(1:npl)%parent = [(i, i=1, npl)] - do i = 1, npl - if (allocated(pl%kin(i)%child)) deallocate(pl%kin(i)%child) - end do + call pl%reset_kinship([(i, i=1, npl)]) ! Re-build the zero-level encounter list, being sure to save the original level information for all bodies allocate(levelg_orig_pl, source=pl%levelg) @@ -540,6 +536,28 @@ module subroutine symba_util_rearray_pl(self, system, param) end subroutine symba_util_rearray_pl + module subroutine symba_util_reset_kinship(self, idx) + !! author: David A. Minton + !! + !! Resets the kinship status of bodies. + !! + implicit none + class(symba_pl), intent(inout) :: self !! SyMBA massive body object + integer(I4B), dimension(:), intent(in) :: idx !! Index array of bodies to reset + ! Internals + integer(I4B) :: i, j + + self%kin(idx(:))%parent = idx(:) + self%kin(idx(:))%nchild = 0 + do j = 1, size(idx(:)) + i = idx(j) + if (allocated(self%kin(i)%child)) deallocate(self%kin(i)%child) + end do + + return + end subroutine symba_util_reset_kinship + + module subroutine symba_util_resize_arr_kin(arr, nnew) !! author: David A. Minton !!