From f378e69719002957f6afbbbb251b84cd215ecabf Mon Sep 17 00:00:00 2001 From: David A Minton Date: Thu, 22 Dec 2022 17:16:11 -0500 Subject: [PATCH] Cleanup and bugfixes --- examples/Fragmentation/Fragmentation_Movie.py | 2 +- src/collision/collision_generate.f90 | 46 ++++++++++++++----- src/collision/collision_module.f90 | 4 +- src/collision/collision_util.f90 | 16 +++++-- 4 files changed, 50 insertions(+), 18 deletions(-) diff --git a/examples/Fragmentation/Fragmentation_Movie.py b/examples/Fragmentation/Fragmentation_Movie.py index baff5c002..b8fd8217a 100644 --- a/examples/Fragmentation/Fragmentation_Movie.py +++ b/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") diff --git a/src/collision/collision_generate.f90 b/src/collision/collision_generate.f90 index d689c3662..421cf0558 100644 --- a/src/collision/collision_generate.f90 +++ b/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 diff --git a/src/collision/collision_module.f90 b/src/collision/collision_module.f90 index 237999adf..7b8331d1b 100644 --- a/src/collision/collision_module.f90 +++ b/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 diff --git a/src/collision/collision_util.f90 b/src/collision/collision_util.f90 index 4474544c1..9338597a3 100644 --- a/src/collision/collision_util.f90 +++ b/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")