Cleanup and bugfixes

examples/Fragmentation/Fragmentation_Movie.py

@@ -206,7 +206,7 @@ def data_stream(self, frame=0):
         # Set fragmentation parameters
         minimum_fragment_gmass = 0.2 * body_Gmass[style][1] # Make the minimum fragment mass a fraction of the smallest body
         gmtiny = 0.99 * body_Gmass[style][1] # Make GMTINY just smaller than the smallest original body. This will prevent runaway collisional cascades
-        sim.set_parameter(collision_model="merge", encounter_save="both", gmtiny=gmtiny, minimum_fragment_gmass=minimum_fragment_gmass, verbose=False)
+        sim.set_parameter(collision_model="bounce", encounter_save="both", gmtiny=gmtiny, minimum_fragment_gmass=minimum_fragment_gmass, verbose=False)
         sim.run(dt=1e-3, tstop=1.0e-3, istep_out=1, dump_cadence=0)
 
         print("Generating animation")

src/collision/collision_generate.f90

@@ -33,10 +33,7 @@ module subroutine collision_generate_hitandrun(collider, nbody_system, param, t)
       class is (swiftest_nbody_system)
       select type (pl => nbody_system%pl)
       class is (swiftest_pl)
-      select type(before => collider%after)
-      class is (swiftest_nbody_system)
-      select type(after => collider%after)
-      class is (swiftest_nbody_system)
+
          associate(impactors => collider%impactors, fragments => collider%fragments,  pl => nbody_system%pl)
             message = "Hit and run between"
             call collision_io_collider_message(nbody_system%pl, impactors%id, message)
@@ -46,12 +43,19 @@ module subroutine collision_generate_hitandrun(collider, nbody_system, param, t)
             pl%status(impactors%id(:)) = ACTIVE
             pl%ldiscard(impactors%id(:)) = .false.
             pl%lcollision(impactors%id(:)) = .false.
-            allocate(after%pl, source=before%pl) ! Be sure to save the pl so that snapshots still work 
+
+            ! Be sure to save the pl so that snapshots still work 
+            select type(before => collider%before)
+            class is (swiftest_nbody_system)
+            select type(after => collider%after)
+            class is (swiftest_nbody_system)
+               allocate(after%pl, source=before%pl) 
+            end select
+            end select
          end associate
       end select
       end select
-      end select
-      end select
+
 
 
       return
@@ -167,17 +171,36 @@ module subroutine collision_generate_bounce(self, nbody_system, param, t)
       class(base_parameters),   intent(inout) :: param        !! Current run configuration parameters 
       real(DP),                 intent(in)    :: t            !! The time of the collision
       ! Internals
-      integer(I4B) :: i,nfrag
+      integer(I4B) :: i,j,nfrag
+      real(DP), dimension(NDIM) :: vcom, rnorm
 
       select type(nbody_system)
       class is (swiftest_nbody_system)
+      select type (pl => nbody_system%pl)
+      class is (swiftest_pl)
          associate(impactors => nbody_system%collider%impactors, fragments => nbody_system%collider%fragments)
             select case (impactors%regime) 
             case (COLLRESOLVE_REGIME_DISRUPTION, COLLRESOLVE_REGIME_SUPERCATASTROPHIC)
                nfrag = size(impactors%id(:))
-               call self%setup_fragments(nfrag)
-
-
+               do i = 1, nfrag
+                  j = impactors%id(i)
+                  vcom(:) = pl%vb(:,j) - impactors%vbcom(:)
+                  rnorm(:) = .unit. (impactors%rb(:,2) - impactors%rb(:,1))
+                  ! Do the reflection
+                  vcom(:) = vcom(:) - 2 * dot_product(vcom(:),rnorm(:)) * rnorm(:)
+                  pl%vb(:,j) = impactors%vbcom(:) + vcom(:)
+                  self%status = DISRUPTED
+                  pl%status(j) = ACTIVE
+                  pl%ldiscard(j) = .false.
+                  pl%lcollision(j) = .false.
+               end do
+               select type(before => self%before)
+               class is (swiftest_nbody_system)
+               select type(after => self%after)
+               class is (swiftest_nbody_system)
+                  allocate(after%pl, source=before%pl) ! Be sure to save the pl so that snapshots still work   
+               end select
+               end select
             case (COLLRESOLVE_REGIME_HIT_AND_RUN)
                call collision_generate_hitandrun(self, nbody_system, param, t)
             case (COLLRESOLVE_REGIME_MERGE, COLLRESOLVE_REGIME_GRAZE_AND_MERGE)
@@ -188,6 +211,7 @@ module subroutine collision_generate_bounce(self, nbody_system, param, t)
             end select
             end associate
       end select
+      end select
 
       return
    end subroutine collision_generate_bounce

src/collision/collision_module.f90

@@ -78,7 +78,7 @@ module collision
       procedure :: consolidate => collision_resolve_consolidate_impactors !! Consolidates a multi-body collision into an equivalent 2-body collision
       procedure :: get_regime  => collision_regime_impactors              !! Determine which fragmentation regime the set of impactors will be
       procedure :: reset       => collision_util_reset_impactors          !! Resets the collider object variables to 0 and deallocates the index and mass distributions
-      final     ::               collision_final_impactors                !! Finalizer will deallocate all allocatables
+      final     ::                collision_final_impactors                !! Finalizer will deallocate all allocatables
    end type collision_impactors
 
 
@@ -353,7 +353,7 @@ end subroutine collision_util_set_coordinate_system
 
       module subroutine collision_setup_system(self, nbody_system)
          implicit none
-         class(collision_merge),  intent(inout) :: self         !! Encounter collision system object
+         class(collision_merge),   intent(inout) :: self         !! Encounter collision system object
          class(base_nbody_system), intent(in)    :: nbody_system !! Current nbody system. Used as a mold for the before/after snapshots
       end subroutine collision_setup_system
 

src/collision/collision_util.f90

@@ -351,8 +351,16 @@ module subroutine collision_util_reset_system(self)
       ! Arguments
       class(collision_merge),    intent(inout) :: self  !! Collision system object
 
-      if (allocated(self%before)) deallocate(self%before)
-      if (allocated(self%after)) deallocate(self%after)
+      select type(before => self%before)
+      class is (swiftest_nbody_system)
+         if (allocated(before%pl)) deallocate(before%pl)
+         if (allocated(before%tp)) deallocate(before%tp)
+      end select
+      select type(after => self%after)
+      class is (swiftest_nbody_system)
+         if (allocated(after%pl)) deallocate(after%pl)
+         if (allocated(after%tp)) deallocate(after%tp)
+      end select
 
       self%Lorbit(:,:) = 0.0_DP
       self%Lspin(:,:) = 0.0_DP
@@ -486,7 +494,7 @@ module subroutine collision_util_snapshot(self, param, nbody_system, t, arg)
       character(*),                intent(in),   optional :: arg    !! "before": takes a snapshot just before the collision. "after" takes the snapshot just after the collision.
       ! Arguments
       class(collision_snapshot), allocatable :: snapshot
-      class(swiftest_pl), allocatable                    :: pl
+      class(swiftest_pl), allocatable :: pl
       character(len=:), allocatable :: stage
 
       if (present(arg)) then
@@ -515,7 +523,7 @@ module subroutine collision_util_snapshot(self, param, nbody_system, t, arg)
                pl%info(:) = nbody_system%pl%info(idx(:))
                select type (before => nbody_system%collider%before)
                class is (swiftest_nbody_system)
-                  allocate(before%pl, source=pl)
+                  call move_alloc(pl, before%pl)
                end select
             end associate
          case("after")