Added merger code

src/modules/rmvs_classes.f90

@@ -70,7 +70,7 @@ module rmvs_classes
       procedure :: encounter_check => rmvs_encounter_check_tp     !! Checks if any test particles are undergoing a close encounter with a massive body
       procedure :: accel           => rmvs_kick_getacch_tp        !! Calculates either the standard or modified version of the acceleration depending if the
                                                                   !!    if the test particle is undergoing a close encounter or not
-      procedure :: setup           => rmvs_setup_tp               !! Constructor method - Allocates space for number of particles
+      procedure :: setup           => rmvs_setup_tp               !! Constructor method - Allocates space for the input number of bodiess
       procedure :: append          => rmvs_util_append_tp         !! Appends elements from one structure to another
       procedure :: fill            => rmvs_util_fill_tp           !! "Fills" bodies from one object into another depending on the results of a mask (uses the UNPACK intrinsic)
       procedure :: resize          => rmvs_util_resize_tp         !! Checks the current size of a Swiftest body against the requested size and resizes it if it is too small.
@@ -93,7 +93,7 @@ module rmvs_classes
       class(rmvs_nbody_system), dimension(:), allocatable :: planetocentric            !! Planetocentric version of the massive body objects (one for each massive body)
       logical                                             :: lplanetocentric = .false. !! Flag that indicates that the object is a planetocentric set of masive bodies used for close encounter calculations
    contains
-      procedure :: setup     => rmvs_setup_pl               !! Constructor method - Allocates space for number of particles
+      procedure :: setup     => rmvs_setup_pl               !! Constructor method - Allocates space for the input number of bodiess
       procedure :: append    => rmvs_util_append_pl         !! Appends elements from one structure to another
       procedure :: fill      => rmvs_util_fill_pl           !! "Fills" bodies from one object into another depending on the results of a mask (uses the UNPACK intrinsic)
       procedure :: resize    => rmvs_util_resize_pl         !! Checks the current size of a Swiftest body against the requested size and resizes it if it is too small.

src/modules/symba_classes.f90

@@ -91,16 +91,24 @@ module symba_classes
       procedure :: drift           => symba_drift_pl                 !! Method for Danby drift in Democratic Heliocentric coordinates. Sets the mask to the current recursion level
       procedure :: encounter_check => symba_encounter_check_pl       !! Checks if massive bodies are going through close encounters with each other
       procedure :: accel           => symba_kick_getacch_pl          !! Compute heliocentric accelerations of massive bodies
-      procedure :: setup           => symba_setup_pl                 !! Constructor method - Allocates space for number of particle
+      procedure :: setup           => symba_setup_pl                 !! Constructor method - Allocates space for the input number of bodies
       procedure :: append          => symba_util_append_pl           !! Appends elements from one structure to another
       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 :: resize          => symba_util_resize_pl           !! Checks the current size of a Swiftest body against the requested size and resizes it if it is too small.
+      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
       procedure :: spill           => symba_util_spill_pl            !! "Spills" bodies from one object to another depending on the results of a mask (uses the PACK intrinsic)
    end type symba_pl
 
+   type, extends(symba_pl) :: symba_merger
+      integer(I4B), dimension(:), allocatable :: ncomp
+   contains
+      procedure :: append          => symba_util_append_merger       !! Appends elements from one structure to another
+      procedure :: resize          => symba_util_resize_merger       !! Checks the current size of a SyMBA merger list against the requested size and resizes it if it is too small.
+      procedure :: setup           => symba_setup_merger             !! Constructor method - Allocates space for the input number of bodies
+   end type symba_merger
+
    !********************************************************************************************************************************
    !                                    symba_tp class definitions and method interfaces
    !*******************************************************************************************************************************
@@ -113,7 +121,7 @@ module symba_classes
       procedure :: drift           => symba_drift_tp               !! Method for Danby drift in Democratic Heliocentric coordinates. Sets the mask to the current recursion level
       procedure :: encounter_check => symba_encounter_check_tp     !! Checks if any test particles are undergoing a close encounter with a massive body
       procedure :: accel           => symba_kick_getacch_tp        !! Compute heliocentric accelerations of test particles
-      procedure :: setup           => symba_setup_tp               !! Constructor method - Allocates space for number of particle
+      procedure :: setup           => symba_setup_tp               !! Constructor method - Allocates space for the input number of bodies
       procedure :: append          => symba_util_append_tp         !! Appends elements from one structure to another
       procedure :: fill            => symba_util_fill_tp           !! "Fills" bodies from one object into another depending on the results of a mask (uses the UNPACK intrinsic)
       procedure :: resize          => symba_util_resize_tp         !! Checks the current size of a Swiftest body against the requested size and resizes it if it is too small.
@@ -151,8 +159,8 @@ module symba_classes
    !  symba_nbody_system class definitions and method interfaces
    !********************************************************************************************************************************
    type, extends(helio_nbody_system) :: symba_nbody_system
-      class(symba_pl),      allocatable :: mergeadd_list !! List of added bodies in mergers or collisions
-      class(symba_pl),      allocatable :: mergesub_list !! List of subtracted bodies in mergers or collisions
+      class(symba_merger),  allocatable :: mergeadd_list !! List of added bodies in mergers or collisions
+      class(symba_merger),  allocatable :: mergesub_list !! List of subtracted bodies in mergers or collisions
       class(symba_pltpenc), allocatable :: pltpenc_list  !! List of massive body-test particle encounters in a single step 
       class(symba_plplenc), allocatable :: plplenc_list  !! List of massive body-massive body encounters in a single step
       integer(I4B)                      :: irec          !! System recursion level
@@ -266,6 +274,16 @@ module function symba_encounter_check_tp(self, system, dt, irec) result(lany_enc
          logical                                  :: lany_encounter !! Returns true if there is at least one close encounter      
       end function symba_encounter_check_tp
 
+      module function symba_fragmentation_casemerge(system, param, family, x, v, mass, radius, L_spin, Ip)  result(status)
+         implicit none
+         class(symba_nbody_system),       intent(inout) :: system           !! SyMBA nbody system object
+         class(symba_parameters),         intent(in)    :: param            !! Current run configuration parameters with SyMBA additions
+         integer(I4B),    dimension(:),   intent(in)    :: family           !! List of indices of all bodies inovlved in the collision
+         real(DP),        dimension(:,:), intent(in)    :: x, v, L_spin, Ip !! Input values that represent a 2-body equivalent of a possibly 2+ body collision
+         real(DP),        dimension(:),   intent(in)    :: mass, radius     !! Input values that represent a 2-body equivalent of a possibly 2+ body collisio
+         integer(I4B)                                   :: status           !! Status flag assigned to this outcome
+      end function symba_fragmentation_casemerge
+
       module subroutine symba_io_write_discard(self, param)
          use swiftest_classes, only : swiftest_parameters
          implicit none
@@ -357,11 +375,26 @@ module subroutine symba_io_write_frame_info(self, iu, param)
          class(swiftest_parameters), intent(in)    :: param !! Current run configuration parameters 
       end subroutine symba_io_write_frame_info 
 
+      module subroutine symba_setup_initialize_system(self, param)
+         use swiftest_classes, only : swiftest_parameters
+         implicit none
+         class(symba_nbody_system),  intent(inout) :: self  !! SyMBA system object
+         class(swiftest_parameters), intent(inout) :: param !! Current run configuration parameters 
+      end subroutine symba_setup_initialize_system
+
+      module subroutine symba_setup_merger(self, n, param)
+         use swiftest_classes, only : swiftest_parameters
+         implicit none
+         class(symba_merger),        intent(inout) :: self  !! SyMBA merger list object
+         integer(I4B),               intent(in)    :: n     !! Number of particles to allocate space for
+         class(swiftest_parameters), intent(in)    :: param !! Current run configuration parameters
+      end subroutine symba_setup_merger
+
       module subroutine symba_setup_pl(self, n, param)
          use swiftest_classes, only : swiftest_parameters
          implicit none
-         class(symba_pl),           intent(inout) :: self  !! SyMBA massive body object
-         integer(I4B),              intent(in)    :: n     !! Number of particles to allocate space for
+         class(symba_pl),            intent(inout) :: self  !! SyMBA massive body object
+         integer(I4B),               intent(in)    :: n     !! Number of particles to allocate space for
          class(swiftest_parameters), intent(in)    :: param !! Current run configuration parameters
       end subroutine symba_setup_pl
 
@@ -371,19 +404,6 @@ module subroutine symba_setup_pltpenc(self,n)
          integer(I4B),         intent(in)    :: n    !! Number of encounters to allocate space for
       end subroutine symba_setup_pltpenc
 
-      module subroutine symba_setup_plplenc(self,n)
-         implicit none
-         class(symba_plplenc), intent(inout) :: self !! SyMBA pl-tp encounter structure
-         integer(I4B),         intent(in)    :: n    !! Number of encounters to allocate space for
-      end subroutine symba_setup_plplenc
-
-      module subroutine symba_setup_initialize_system(self, param)
-         use swiftest_classes, only : swiftest_parameters
-         implicit none
-         class(symba_nbody_system),  intent(inout) :: self  !! SyMBA system object
-         class(swiftest_parameters), intent(inout) :: param !! Current run configuration parameters 
-      end subroutine symba_setup_initialize_system
-
       module subroutine symba_setup_tp(self, n, param)
          use swiftest_classes, only : swiftest_parameters
          implicit none
@@ -448,6 +468,14 @@ end subroutine symba_util_append_arr_kin
    end interface
 
    interface
+      module subroutine symba_util_append_merger(self, source, lsource_mask)
+         use swiftest_classes, only : swiftest_body
+         implicit none
+         class(symba_merger),             intent(inout) :: self         !! SyMBA massive body object
+         class(swiftest_body),            intent(in)    :: source       !! Source object to append
+         logical, dimension(:), optional, intent(in)    :: lsource_mask !! Logical mask indicating which elements to append to
+      end subroutine symba_util_append_merger
+
       module subroutine symba_util_append_pl(self, source, lsource_mask)
          use swiftest_classes, only : swiftest_body
          implicit none
@@ -522,6 +550,12 @@ end subroutine symba_util_resize_arr_kin
    end interface
 
    interface
+      module subroutine symba_util_resize_merger(self, nnew)
+         implicit none
+         class(symba_merger), intent(inout) :: self  !! SyMBA merger list object
+         integer(I4B),        intent(in)    :: nnew  !! New size neded
+      end subroutine symba_util_resize_merger
+
       module subroutine symba_util_resize_pl(self, nnew)
          implicit none
          class(symba_pl), intent(inout) :: self  !! SyMBA massive body object

src/modules/whm_classes.f90

@@ -46,7 +46,7 @@ module whm_classes
       procedure :: sort        => whm_util_sort_pl           !! Sort a WHM massive body object in-place. 
       procedure :: rearrange   => whm_util_sort_rearrange_pl !! Rearranges the order of array elements of body based on an input index array. Used in sorting methods
       procedure :: spill       => whm_util_spill_pl          !!"Spills" bodies from one object to another depending on the results of a mask (uses the PACK intrinsic)
-      procedure :: setup       => whm_setup_pl               !! Constructor method - Allocates space for number of particles
+      procedure :: setup       => whm_setup_pl               !! Constructor method - Allocates space for the input number of bodiess
       procedure :: step        => whm_step_pl                !! Steps the body forward one stepsize
    end type whm_pl
 

src/setup/setup.f90

@@ -54,8 +54,8 @@ module subroutine setup_construct_system(system, param)
             allocate(symba_pl :: system%pl)
             allocate(symba_tp :: system%tp)
             allocate(symba_tp :: system%tp_discards)
-            allocate(symba_pl :: system%mergeadd_list)
-            allocate(symba_pl :: system%mergesub_list)
+            allocate(symba_merger :: system%mergeadd_list)
+            allocate(symba_merger :: system%mergesub_list)
             allocate(symba_plplenc :: system%plplenc_list)
             allocate(symba_pltpenc :: system%pltpenc_list)
          end select

src/symba/symba_collision.f90

@@ -208,8 +208,8 @@ function symba_collision_consolidate_familes(pl, param, idx_parent, family, x, v
       fam_size = 2 + sum(nchild(:))
       allocate(family(fam_size))
       family = [parent_child_index_array(1)%idx(:),parent_child_index_array(2)%idx(:)]
-      fam_size = count(pl%status(family(:)) == ACTIVE)
-      family = pack(family(:), pl%status(family(:)) == ACTIVE)
+      fam_size = count(pl%status(family(:)) == COLLISION)
+      family = pack(family(:), pl%status(family(:)) == COLLISION)
       L_spin(:,:) = 0.0_DP
       Ip(:,:) = 0.0_DP
 
@@ -423,12 +423,12 @@ module subroutine symba_collision_resolve_mergers(self, system, param)
       class(symba_nbody_system), intent(inout) :: system !! SyMBA nbody system object
       class(symba_parameters),   intent(in)    :: param  !! Current run configuration parameters with SyMBA additions
       ! Internals
-      integer(I4B) :: i
-      logical      :: lgoodcollision
-      integer(I4B),    dimension(:),     allocatable :: family           !! List of indices of all bodies inovlved in the collision
-      integer(I4B),    dimension(2)                  :: idx_parent       !! Index of the two bodies considered the "parents" of the collision
-      real(DP),        dimension(NDIM,2)             :: x, v, L_spin, Ip !! Output values that represent a 2-body equivalent of a possibly 2+ body collision
-      real(DP),        dimension(2)                  :: mass, radius     !! Output values that represent a 2-body equivalent of a possibly 2+ body collision
+      integer(I4B), dimension(:),     allocatable :: family           !! List of indices of all bodies inovlved in the collision
+      integer(I4B), dimension(2)                  :: idx_parent       !! Index of the two bodies considered the "parents" of the collision
+      real(DP),     dimension(NDIM,2)             :: x, v, L_spin, Ip !! Output values that represent a 2-body equivalent of a possibly 2+ body collision
+      real(DP),     dimension(2)                  :: mass, radius     !! Output values that represent a 2-body equivalent of a possibly 2+ body collision
+      logical                                     :: lgoodcollision
+      integer(I4B)                                :: i, status
 
       associate(plpl_collisions => self, ncollisions => self%nenc, idx1 => self%index1, idx2 => self%index2)
          select type(pl => system%pl)
@@ -439,7 +439,7 @@ module subroutine symba_collision_resolve_mergers(self, system, param)
                lgoodcollision = symba_collision_consolidate_familes(pl, param, idx_parent, family, x, v, mass, radius, L_spin, Ip)
                if (.not. lgoodcollision) cycle
                if (any(pl%status(idx_parent(:)) /= COLLISION)) cycle ! One of these two bodies has already been resolved
-               !call symba_collision_casemerge(system, param, family, x, v, mass, radius, L_spin, Ip) 
+               status = symba_fragmentation_casemerge(system, param, family, x, v, mass, radius, L_spin, Ip) 
             end do
          end select
       end associate

src/symba/symba_discard.f90

@@ -153,6 +153,7 @@ module subroutine symba_discard_pl(self, system, param)
                call symba_discard_nonplpl(self, system, param)
                call plplenc_list%scrub_non_collision(system, param)
                if (plplenc_list%nenc == 0) return ! No collisions to resolve
+               write(*, *) "Collision detected at time t = ",param%t
 
                call pl%h2b(system%cb) 
                if (param%lfragmentation) then