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

Commit

Permalink
Added encounter checks for the plpl and pltp encounter lists
Browse files Browse the repository at this point in the history
  • Loading branch information
@daminton
daminton committed Jul 24, 2021
1 parent d18ffb5 commit 43f6b4c
Show file tree
Hide file tree
Showing 3 changed files with 119 additions and 11 deletions.
20 changes: 20 additions & 0 deletions src/modules/symba_classes.f90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -121,6 +121,7 @@ module symba_classes
integer(I4B), dimension(:), allocatable :: index1 !! position of the planet in encounter
integer(I4B), dimension(:), allocatable :: index2 !! position of the test particle in encounter
contains
procedure, public :: encounter_check => symba_encounter_check_pltpenc !! Checks if massive bodies are going through close encounters with each other
procedure, public :: setup => symba_setup_pltpenc !! A constructor that sets the number of encounters and allocates and initializes all arrays
procedure, public :: copy => symba_util_copy_pltpenc !! Copies all elements of one pltpenc list to another
procedure, public :: resize => symba_util_resize_pltpenc !! Checks the current size of the pltpenc_list against the required size and extends it by a factor of 2 more than requested if it is too small
Expand All @@ -136,6 +137,7 @@ module symba_classes
real(DP), dimension(:,:), allocatable :: vb1 !! the barycentric velocity of parent 1 in encounter
real(DP), dimension(:,:), allocatable :: vb2 !! the barycentric velocity of parent 2 in encounter
contains
procedure, public :: encounter_check => symba_encounter_check_plplenc !! Checks if massive bodies are going through close encounters with each other
procedure, public :: setup => symba_setup_plplenc !! A constructor that sets the number of encounters and allocates and initializes all arrays
procedure, public :: copy => symba_util_copy_plplenc !! Copies all elements of one plplenc list to another
end type symba_plplenc
Expand Down Expand Up @@ -192,6 +194,24 @@ module function symba_encounter_check_pl(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_pl

module function symba_encounter_check_plplenc(self, system, dt, irec) result(lany_encounter)
implicit none
class(symba_plplenc), intent(inout) :: self !! SyMBA pl-pl encounter list object
class(symba_nbody_system), intent(inout) :: system !! SyMBA nbody system object
real(DP), intent(in) :: dt !! step size
integer(I4B), intent(in) :: irec !! Current recursion level
logical :: lany_encounter !! Returns true if there is at least one close encounter
end function symba_encounter_check_plplenc

module function symba_encounter_check_pltpenc(self, system, dt, irec) result(lany_encounter)
implicit none
class(symba_pltpenc), intent(inout) :: self !! SyMBA pl-pl encounter list object
class(symba_nbody_system), intent(inout) :: system !! SyMBA nbody system object
real(DP), intent(in) :: dt !! step size
integer(I4B), intent(in) :: irec !! Current recursion level
logical :: lany_encounter !! Returns true if there is at least one close encounter
end function symba_encounter_check_pltpenc

module function symba_encounter_check_tp(self, system, dt, irec) result(lany_encounter)
implicit none
class(symba_tp), intent(inout) :: self !! SyMBA test particle object
Expand Down
96 changes: 96 additions & 0 deletions src/symba/symba_encounter_check.f90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -50,6 +50,102 @@ module function symba_encounter_check_pl(self, system, dt, irec) result(lany_enc
return
end function symba_encounter_check_pl

module function symba_encounter_check_plplenc(self, system, dt, irec) result(lany_encounter)
!! author: David A. Minton
!!
!! Check for an encounter between massive bodies in the plplenc list.
!!
!! Adapted from portions of David E. Kaufmann's Swifter routine: symba_step_recur.f90
implicit none
! Arguments
class(symba_plplenc), intent(inout) :: self !! SyMBA pl-pl encounter list object
class(symba_nbody_system), intent(inout) :: system !! SyMBA nbody system object
real(DP), intent(in) :: dt !! step size
integer(I4B), intent(in) :: irec !! Current recursion level
logical :: lany_encounter !! Returns true if there is at least one close encounter
! Internals
integer(I4B) :: i
real(DP), dimension(NDIM) :: xr, vr
logical :: lencounter
real(DP) :: rlim2, rji2

lany_encounter = .false.
associate(plplenc_list => self, nplplenc => self%nenc)
select type(pl => system%pl)
class is (symba_pl)
do i = 1, nplplenc
associate(index_i => plplenc_list%index1(i), index_j => plplenc_list%index2(i))
xr(:) = pl%xh(:,index_j) - pl%xh(:,index_i)
vr(:) = pl%vb(:,index_j) - pl%vb(:,index_i)
call symba_encounter_check_one(xr(1), xr(2), xr(3), vr(1), vr(2), vr(3), pl%rhill(index_i), pl%rhill(index_j), dt, irec, lencounter, plplenc_list%lvdotr(i))
if (lencounter) then
rlim2 = (pl%radius(index_i) + pl%radius(index_j))**2
rji2 = dot_product(xr(:), xr(:))! Check to see if these are physically overlapping bodies first, which we should ignore
if (rji2 > rlim2) then
lany_encounter = .true.
pl%levelg(index_i) = irec
pl%levelm(index_i) = MAX(irec, pl%levelm(index_i))
pl%levelg(index_j) = irec
pl%levelm(index_j) = MAX(irec, pl%levelm(index_j))
plplenc_list%level(i) = irec
end if
end if
end associate
end do
end select
end associate
return
end function symba_encounter_check_plplenc

module function symba_encounter_check_pltpenc(self, system, dt, irec) result(lany_encounter)
!! author: David A. Minton
!!
!! Check for an encounter between test particles and massive bodies in the pltpenc list.
!!
!! Adapted from portions of David E. Kaufmann's Swifter routine: symba_step_recur.f90
implicit none
! Arguments
class(symba_pltpenc), intent(inout) :: self !! SyMBA pl-pl encounter list object
class(symba_nbody_system), intent(inout) :: system !! SyMBA nbody system object
real(DP), intent(in) :: dt !! step size
integer(I4B), intent(in) :: irec !! Current recursion level
logical :: lany_encounter !! Returns true if there is at least one close encounter
! Internals
integer(I4B) :: i
real(DP), dimension(NDIM) :: xr, vr
logical :: lencounter
real(DP) :: rlim2, rji2

lany_encounter = .false.
associate(pltpenc_list => self, npltpenc => self%nenc)
select type(pl => system%pl)
class is (symba_pl)
select type(tp => system%tp)
class is (symba_tp)
do i = 1, npltpenc
associate(index_i => pltpenc_list%index1(i), index_j => pltpenc_list%index2(i))
xr(:) = tp%xh(:,index_j) - pl%xh(:,index_i)
vr(:) = tp%vb(:,index_j) - pl%vb(:,index_i)
call symba_encounter_check_one(xr(1), xr(2), xr(3), vr(1), vr(2), vr(3), pl%rhill(index_i), 0.0_DP, dt, irec, lencounter, pltpenc_list%lvdotr(i))
if (lencounter) then
rlim2 = (pl%radius(index_i))**2
rji2 = dot_product(xr(:), xr(:))! Check to see if these are physically overlapping bodies first, which we should ignore
if (rji2 > rlim2) then
lany_encounter = .true.
pl%levelg(index_i) = irec
pl%levelm(index_i) = MAX(irec, pl%levelm(index_i))
tp%levelg(index_j) = irec
tp%levelm(index_j) = MAX(irec, tp%levelm(index_j))
pltpenc_list%level(i) = irec
end if
end if
end associate
end do
end select
end select
end associate
end function symba_encounter_check_pltpenc

module function symba_encounter_check_tp(self, system, dt, irec) result(lany_encounter)
!! author: David A. Minton
!!
Expand Down
14 changes: 3 additions & 11 deletions src/symba/symba_step.f90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -112,12 +112,12 @@ module recursive subroutine symba_step_recur_system(self, param, t, dt, ireci)
real(DP), intent(in) :: dt !! Current stepsize
integer(I4B), value, intent(in) :: ireci !! input recursion level
! Internals
integer(I4B) :: i, j, irecp, icflg, nloops
integer(I4B) :: i, j, irecp, nloops
real(DP) :: dtl, dth, sgn
real(DP), dimension(NDIM) :: xr, vr
logical :: lencounter

associate(pl => self%pl, tp => self%tp, plplenc_list => self%plplenc_list, nplplenc => self%plplenc_list%nenc, pltpenc_list => self%pltpenc_list, npltpenc => self%pltpenc_list%nenc)
associate(system => self, pl => self%pl, tp => self%tp, plplenc_list => self%plplenc_list, pltpenc_list => self%pltpenc_list)
dtl = param%dt / (NTENC**ireci)
dth = 0.5_DP * dtl
IF (dtl / param%dt < VSMALL) THEN
Expand All @@ -133,15 +133,7 @@ module recursive subroutine symba_step_recur_system(self, param, t, dt, ireci)
nloops = NTENC
end if
do j = 1, nloops
icflg = 0
do i = 1, nplplenc
associate(index_i => plplenc_list%index1(i), index_j => plplenc_list%index2(i))
xr(:) = pl%xh(:,index_j) - pl%xh(:,index_i)
vr(:) = pl%vb(:,index_j) - pl%vb(:,index_i)
call symba_encounter_check_one(xr(1), xr(2), xr(3), vr(1), vr(2), vr(3), pl%rhill(index_i), pl%rhill(index_j), dtl, irecp, lencounter, plplenc_list%lvdotr(i))
end associate
end do

lencounter = plplenc_list%encounter_check(system, dtl, irecp) .or. pltpenc_list%encounter_check(system, dtl, irecp)
end do
end associate

Expand Down

0 comments on commit 43f6b4c

Please sign in to comment.