Skip to content
This repository was archived by the owner on Aug 28, 2024. It is now read-only.

Commit

Permalink
Refactored to enforce line limit and started started to switch to tra…
Browse files Browse the repository at this point in the history 
…cking all discards through the collision module
  • Loading branch information
@daminton
daminton committed Feb 22, 2024
1 parent d96884a commit a6a96e9
Show file tree
Hide file tree
Showing 3 changed files with 68 additions and 42 deletions.
62 changes: 36 additions & 26 deletions src/collision/collision_resolve.f90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -11,19 +11,19 @@
use swiftest
contains


module subroutine collision_resolve_consolidate_impactors(self, nbody_system, param, idx_parent, lflag)
!! author: David A. Minton
!!
!! Loops through the pl-pl collision list and groups families together by index. Outputs the indices of all impactors%id members,
!! and pairs of quantities (x and v vectors, mass, radius, L_spin, and Ip) that can be used to resolve the collisional outcome.
!! Loops through the pl-pl collision list and groups families together by index. Outputs the indices of all impactors%id
!! members, and pairs of quantities (x and v vectors, mass, radius, L_spin, and Ip) that can be used to resolve the
!! collisional outcome.
implicit none
! Arguments
class(collision_impactors), intent(out) :: self !! Collision impactors object
class(base_nbody_system), intent(inout) :: nbody_system !! Swiftest nbody system object
class(base_parameters), intent(in) :: param !! Current run configuration parameters with Swiftest additions
integer(I4B), dimension(:), intent(inout) :: idx_parent !! Index of the two bodies considered the "parents" of the collision
logical, intent(out) :: lflag !! Logical flag indicating whether a impactors%id was successfully created or not
class(collision_impactors),intent(out) :: self !! Collision impactors object
class(base_nbody_system),intent(inout) :: nbody_system !! Swiftest nbody system object
class(base_parameters), intent(in) :: param !! Current run configuration parameters with Swiftest additions
integer(I4B), dimension(:), intent(inout) :: idx_parent !! Index of the two bodies considered the "parents" of the collision
logical, intent(out) :: lflag !! Logical flag indicating whether a impactors%id was successfully created or not
! Internals
type collidx_array
integer(I4B), dimension(:), allocatable :: id
Expand All @@ -45,7 +45,8 @@ module subroutine collision_resolve_consolidate_impactors(self, nbody_system, pa
! If all of these bodies share a parent, but this is still a unique collision, move the last child
! out of the parent's position and make it the secondary body
if (idx_parent(1) == idx_parent(2)) then
if (nchild(1) == 0) then ! There is only one valid body recorded in this pair (this could happen due to restructuring of the kinship relationships, though it should be rare)
if (nchild(1) == 0) then ! There is only one valid body recorded in this pair (this could happen due to restructuring
! of the kinship relationships, though it should be rare)
lflag = .false.
call pl%reset_kinship([idx_parent(1)])
return
Expand Down Expand Up @@ -233,15 +234,16 @@ module subroutine collision_resolve_extract_plpl(self, nbody_system, param)
end associate
end do

! Reassemble collision index list to include only those containing the unique pairs of parents, plus all the non-unique pairs that don't
! contain a parent body on the unique parent list.
! Reassemble collision index list to include only those containing the unique pairs of parents, plus all the non-unique
! pairs that don't contain a parent body on the unique parent list.
ncollisions = nunique_parent + count(lplpl_unique_parent)
collision_idx = [unique_parent_idx(:), pack(collision_idx(:), lplpl_unique_parent(:))]

! Create a mask that contains only the pl-pl encounters that did not result in a collision, and then discard them
lplpl_collision(:) = .false.
lplpl_collision(collision_idx(:)) = .true.
call self%spill(nbody_system%plpl_collision, lplpl_collision, ldestructive=.true.) ! Extract any encounters that are not collisions from the list.
! Extract any encounters that are not collisions from the list.
call self%spill(nbody_system%plpl_collision, lplpl_collision, ldestructive=.true.)
end associate
end select
end select
Expand Down Expand Up @@ -281,7 +283,8 @@ module subroutine collision_resolve_extract_pltp(self, nbody_system, param)
! Create a mask that contains only the pl-tp encounters that did not result in a collision, and then discard them
lpltp_collision(:) = .false.
lpltp_collision(collision_idx(:)) = .true.
call self%spill(nbody_system%pltp_collision, lpltp_collision, ldestructive=.true.) ! Extract any encounters that are not collisions from the list.
! Extract any encounters that are not collisions from the list.
call self%spill(nbody_system%pltp_collision, lpltp_collision, ldestructive=.true.)
end associate
end select
end select
Expand Down Expand Up @@ -383,8 +386,9 @@ module subroutine collision_resolve_mergeaddsub(nbody_system, param, t, status)
class is (swiftest_nbody_system)
select type(param)
class is (swiftest_parameters)
associate(pl => nbody_system%pl, pl_discards => nbody_system%pl_discards, info => nbody_system%pl%info, pl_adds => nbody_system%pl_adds, cb => nbody_system%cb, &
collider => nbody_system%collider, impactors => nbody_system%collider%impactors,fragments => nbody_system%collider%fragments)
associate(pl => nbody_system%pl, pl_discards => nbody_system%pl_discards, info => nbody_system%pl%info, &
pl_adds => nbody_system%pl_adds, cb => nbody_system%cb, collider => nbody_system%collider, &
impactors => nbody_system%collider%impactors,fragments => nbody_system%collider%fragments)

! Add the impactors%id bodies to the subtraction list
nimpactors = impactors%ncoll
Expand Down Expand Up @@ -646,10 +650,12 @@ module subroutine collision_resolve_plpl(self, nbody_system, param, t, dt, irec)
allocate(idnew, source=nbody_system%pl_adds%id)
mnew = sum(nbody_system%pl_adds%mass(:))

! Rearrange the arrays: Remove discarded bodies, add any new bodies, re-sort, and recompute all indices and encounter lists
! Rearrange the arrays: Remove discarded bodies, add any new bodies, re-sort, and recompute all indices and
! encounter lists
call pl%rearray(nbody_system, param)

! Destroy the add/discard list so that we don't append the same body multiple times if another collision is detected
! Destroy the add/discard list so that we don't append the same body multiple times if another collision
! is detected
call nbody_system%pl_discards%setup(0, param)
call nbody_system%pl_adds%setup(0, param)

Expand All @@ -673,13 +679,17 @@ module subroutine collision_resolve_plpl(self, nbody_system, param, t, dt, irec)
end if


! Check whether or not any of the particles that were just added are themselves in a collision state. This will generate a new plpl_collision
! Check whether or not any of the particles that were just added are themselves in a collision state. This will
! generate a new plpl_collision
call self%collision_check(nbody_system, param, t, dt, irec, lplpl_collision)

if (.not.lplpl_collision) exit
if (loop == MAXCASCADE) then
call swiftest_io_log_one_message(COLLISION_LOG_OUT,"A runaway collisional cascade has been detected in collision_resolve_plpl.")
call swiftest_io_log_one_message(COLLISION_LOG_OUT,"Consider reducing the step size or changing the parameters in the collisional model to reduce the number of fragments.")
call swiftest_io_log_one_message(COLLISION_LOG_OUT,"A runaway collisional cascade has been detected in " // &
"collision_resolve_plpl.")
call swiftest_io_log_one_message(COLLISION_LOG_OUT,"Consider reducing the step size or changing the " // &
"parameters in the collisional model to reduce the " // &
"number of fragments.")
call base_util_exit(FAILURE,unit=param%display_unit)
end if
end associate
Expand Down Expand Up @@ -740,12 +750,12 @@ module subroutine collision_resolve_pltp(self, nbody_system, param, t, dt, irec)
ncollisions = pltp_collision%nenc
write(timestr,*) t
call swiftest_io_log_one_message(COLLISION_LOG_OUT, "")
call swiftest_io_log_one_message(COLLISION_LOG_OUT, "***********************************************************" // &
"***********************************************************")
call swiftest_io_log_one_message(COLLISION_LOG_OUT, "Collision between test particle and massive body detected at time t = " // &
trim(adjustl(timestr)))
call swiftest_io_log_one_message(COLLISION_LOG_OUT, "***********************************************************" // &
"***********************************************************")
call swiftest_io_log_one_message(COLLISION_LOG_OUT,"***********************************************************" // &
"***********************************************************")
call swiftest_io_log_one_message(COLLISION_LOG_OUT, "Collision between test particle and massive body detected " // &
"at time t = " // trim(adjustl(timestr)))
call swiftest_io_log_one_message(COLLISION_LOG_OUT,"***********************************************************" // &
"***********************************************************")

do k = 1_I8B, ncollisions
! Advance the collision id number and save it
Expand Down
4 changes: 3 additions & 1 deletion src/collision/collision_util.f90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -285,6 +285,7 @@ module subroutine collision_util_dealloc_snapshot(self)
return
end subroutine collision_util_dealloc_snapshot


module subroutine collision_util_dealloc_impactors(self)
!! author: David A. Minton
!!
Expand Down Expand Up @@ -419,7 +420,8 @@ end subroutine collision_util_dealloc_basic
module subroutine collision_util_reset_fragments(self)
!! author: David A. Minton
!!
!! Resets all position and velocity-dependent fragment quantities in order to do a fresh calculation (does not reset mass, radius, or other values that get set prior to the call to fraggle_generate)
!! Resets all position and velocity-dependent fragment quantities in order to do a fresh calculation (does not reset mass,
!! radius, or other values that get set prior to the call to fraggle_generate)
implicit none
! Arguments
class(collision_fragments), intent(inout) :: self
Expand Down
44 changes: 29 additions & 15 deletions src/swiftest/swiftest_discard.f90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -25,7 +25,7 @@ module subroutine swiftest_discard_system(self, param)
lpl_check = allocated(self%pl_discards)
ltp_check = allocated(self%tp_discards)

associate(nbody_system => self, tp => self%tp, pl => self%pl, tp_discards => self%tp_discards, pl_discards => self%pl_discards)
associate(nbody_system => self,tp => self%tp,pl => self%pl,tp_discards => self%tp_discards,pl_discards => self%pl_discards)
lpl_discards = .false.
ltp_discards = .false.
if (lpl_check .and. pl%nbody > 0) then
Expand Down Expand Up @@ -59,8 +59,8 @@ end subroutine swiftest_discard_system
module subroutine swiftest_discard_pl(self, nbody_system, param)
!! author: David A. Minton
!!
!! Placeholder method for discarding massive bodies. This method does nothing except to ensure that the discard flag is set to false.
!! This method is intended to be overridden by more advanced integrators.
!! Placeholder method for discarding massive bodies. This method does nothing except to ensure that the discard flag is set
!! to false. This method is intended to be overridden by more advanced integrators.
implicit none
! Arguments
class(swiftest_pl), intent(inout) :: self !! Swiftest massive body object
Expand All @@ -87,27 +87,40 @@ module subroutine swiftest_discard_tp(self, nbody_system, param)
class(swiftest_nbody_system), intent(inout) :: nbody_system !! Swiftest nbody system object
class(swiftest_parameters), intent(inout) :: param !! Current run configuration parameter
! Internals
logical, dimension(:), allocatable :: ldiscard
integer(I4B) :: npl, ntp
logical, dimension(self%nbody) :: ldiscard
integer(I4B) :: i, nstart, nend, nsub
class(swiftest_tp), allocatable :: tpsub

if (self%nbody == 0) return

associate(tp => self, cb => nbody_system%cb, pl => nbody_system%pl)
ntp = tp%nbody
npl = pl%nbody
associate(tp => self, ntp => self%nbody, cb => nbody_system%cb, pl => nbody_system%pl, npl => nbody_system%pl%nbody, &
tp_discards => nbody_system%tp_discards)

if ((param%rmin >= 0.0_DP) .or. (param%rmax >= 0.0_DP) .or. &
(param%rmaxu >= 0.0_DP) .or. ((param%qmin >= 0.0_DP) .and. (param%qmin_coord == "BARY"))) then
call pl%h2b(cb)
call tp%h2b(cb)
end if

if ((param%rmin >= 0.0_DP) .or. (param%rmax >= 0.0_DP) .or. (param%rmaxu >= 0.0_DP)) call swiftest_discard_cb_tp(tp, nbody_system, param)
if (param%qmin >= 0.0_DP) call swiftest_discard_peri_tp(tp, nbody_system, param)
if (param%lclose) call swiftest_discard_pl_tp(tp, nbody_system, param)
if ((param%rmin >= 0.0_DP) .or. &
(param%rmax >= 0.0_DP) .or. &
(param%rmaxu >= 0.0_DP)) then
call swiftest_discard_cb_tp(tp, nbody_system, param)
end if
if (param%qmin >= 0.0_DP) then
call swiftest_discard_peri_tp(tp, nbody_system, param)
end if
if (param%lclose) then
call swiftest_discard_pl_tp(tp, nbody_system, param)
end if
if (any(tp%ldiscard(1:ntp))) then
allocate(ldiscard, source=tp%ldiscard)
call tp%spill(nbody_system%tp_discards, ldiscard(1:ntp), ldestructive=.true.)
ldiscard(1:ntp) = tp%ldiscard(1:ntp)
allocate(tpsub, mold=tp)
call tp%spill(tpsub, ldiscard, ldestructive=.false.)
nsub = tpsub%nbody
nstart = tp_discards%nbody + 1
nend = tp_discards%nbody + nsub
call tp_discards%append(tpsub, lsource_mask=[(.true., i = 1, nsub)])
end if
end associate

Expand Down Expand Up @@ -275,8 +288,9 @@ subroutine swiftest_discard_pl_tp(tp, nbody_system, param)
write(idstrj, *) pl%id(j)
write(timestr, *) nbody_system%t
write(message, *) "Particle " // trim(adjustl(tp%info(i)%name)) // " (" // trim(adjustl(idstri)) // ")" &
// " too close to massive body " // trim(adjustl(pl%info(j)%name)) // " (" // trim(adjustl(idstrj)) // ")" &
// " at t = " // trim(adjustl(timestr))
// " too close to massive body " // trim(adjustl(pl%info(j)%name)) // " (" &
// trim(adjustl(idstrj)) // ")" &
// " at t = " // trim(adjustl(timestr))
call swiftest_io_log_one_message(COLLISION_LOG_OUT, message)
tp%ldiscard(i) = .true.
call tp%info(i)%set_value(status="DISCARDED_PLR", discard_time=nbody_system%t, discard_rh=tp%rh(:,i), &
Expand Down

0 comments on commit a6a96e9

Please sign in to comment.