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

Commit

Permalink
Browse files Browse the repository at this point in the history 
Simplified RMVS data structures by removing redundant leftovers from the Swifter conversion
  • Loading branch information
@daminton
daminton committed Apr 2, 2021
1 parent d890867 commit e59535c
Show file tree
Hide file tree
Showing 6 changed files with 37 additions and 75 deletions.
3 changes: 1 addition & 2 deletions src/modules/rmvs_classes.f90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -254,7 +254,7 @@ module function rmvs_encounter_check_tp(self, cb, pl, dt, rts) result(lencounter

end function rmvs_encounter_check_tp

module subroutine rmvs_peri_tp(self, cb, pl, t, dt, lfirst, index, nenc, ipleP, config)
module subroutine rmvs_peri_tp(self, cb, pl, t, dt, lfirst, index, ipleP, config)
use swiftest_classes, only : swiftest_configuration
class(rmvs_tp), intent(inout) :: self !! RMVS test particle object
class(rmvs_cb), intent(inout) :: cb !! RMVS central body object
Expand All @@ -263,7 +263,6 @@ module subroutine rmvs_peri_tp(self, cb, pl, t, dt, lfirst, index, nenc, ipleP,
real(DP), intent(in) :: dt !! step size
logical, intent(in) :: lfirst !! Logical flag indicating whether current invocation is the first
integer(I4B), intent(in) :: index !! outer substep number within current set
integer(I4B), intent(in) :: nenc !! number of test particles encountering current planet
integer(I4B), intent(in) :: ipleP !!index of RMVS planet being closely encountered
class(swiftest_configuration), intent(in) :: config !! Input collection of configuration parameters
end subroutine rmvs_peri_tp
Expand Down
17 changes: 0 additions & 17 deletions src/rmvs/rmvs_encounter_check.f90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -30,16 +30,9 @@ module function rmvs_encounter_check_tp(self, cb, pl, dt, rts) result(lencounter
pl%encmask(:,:) = .false.
lencounter = .false.
pl%nenc(:) = 0
pl%tpenc1P(:) = 0
! if first time through, calc hill's sphere for the planets
!if (lfirst) then

! lfirst = .false.
!end if
do i = 1, ntp
if (tp%status(i) == ACTIVE) then
tp%plencP(i) = 0
!tp%tpencP(i) = 0
lflag = .false.

do j = 1, npl
Expand All @@ -51,16 +44,6 @@ module function rmvs_encounter_check_tp(self, cb, pl, dt, rts) result(lencounter
lencounter = .true.
pl%encmask(i,j) = .true.
pl%nenc(j) = pl%nenc(j) + 1
!nenc = pl%nenc(j)
!if (nenc == 1) then
! pl%tpenc1P(j) = i
!else
! tpencPindex = pl%tpenc1P(j)
! do k = 2, nenc - 1
! tpencPindex = tp%tpencP(tpencPindex)
! end do
! tp%tpencP(tpencPindex) = i
!end if
tp%plencP(i) = j
exit
end if
Expand Down
3 changes: 1 addition & 2 deletions src/rmvs/rmvs_peri.f90
View file
Original file line number Diff line number Diff line change
@@ -1,7 +1,7 @@
submodule(rmvs_classes) s_rmvs_peri
use swiftest
contains
module subroutine rmvs_peri_tp(self, cb, pl, t, dt, lfirst, index, nenc, ipleP, config)
module subroutine rmvs_peri_tp(self, cb, pl, t, dt, lfirst, index, ipleP, config)
!! author: David A. Minton
!!
!! Determine planetocentric pericenter passages for test particles in close encounters with a planet
Expand All @@ -17,7 +17,6 @@ module subroutine rmvs_peri_tp(self, cb, pl, t, dt, lfirst, index, nenc, ipleP,
real(DP), intent(in) :: dt !! step size
logical, intent(in) :: lfirst !! Logical flag indicating whether current invocation is the first
integer(I4B), intent(in) :: index !! outer substep number within current set
integer(I4B), intent(in) :: nenc !! number of test particles encountering current planet
integer(I4B), intent(in) :: ipleP !!index of RMVS planet being closely encountered
class(swiftest_configuration), intent(in) :: config !! Input collection of configuration parameters
! Internals
Expand Down
4 changes: 0 additions & 4 deletions src/rmvs/rmvs_setup.f90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -19,7 +19,6 @@ module subroutine rmvs_setup_pl(self,n)
if (n <= 0) return

allocate(self%nenc(n))
allocate(self%tpenc1P(n))
allocate(self%xout(NDIM, n, 0:NTENC))
allocate(self%vout(NDIM, n, 0:NTENC))
allocate(self%xin(NDIM, n, 0:NTPHENC))
Expand All @@ -32,7 +31,6 @@ module subroutine rmvs_setup_pl(self,n)
self%plenc(:)%nbody = n

self%nenc = 0
self%tpenc1P(:) = 0
self%xout(:,:,:) = 0.0_DP
self%vout(:,:,:) = 0.0_DP
self%xin(:,:,:) = 0.0_DP
Expand Down Expand Up @@ -65,12 +63,10 @@ module subroutine rmvs_setup_tp(self,n)
allocate(self%lperi(n))
allocate(self%plperP(n))
allocate(self%plencP(n))
allocate(self%tpencP(n))

self%lperi(:) = .false.
self%plperP(:) = 0
self%plencP(:) = 0
self%tpencP(:) = 0

return
end subroutine rmvs_setup_tp
Expand Down
10 changes: 0 additions & 10 deletions src/rmvs/rmvs_spill_and_fill.f90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -20,10 +20,8 @@ module subroutine rmvs_spill_pl(self, discards, lspill_list)
class is (rmvs_pl)

discards%nenc(:) = pack(keeps%nenc(:), lspill_list(:))
discards%tpenc1P(:) = pack(keeps%tpenc1P(:), lspill_list(:))
if (count(.not.lspill_list(:)) > 0) then
keeps%nenc(:) = pack(keeps%nenc(:), .not. lspill_list(:))
keeps%tpenc1P(:) = pack(keeps%tpenc1P(:), .not. lspill_list(:))
end if
call whm_spill_pl(keeps, discards, lspill_list)
class default
Expand Down Expand Up @@ -56,9 +54,6 @@ module subroutine rmvs_fill_pl(self, inserts, lfill_list)
keeps%nenc(:) = unpack(keeps%nenc(:), .not.lfill_list(:), keeps%nenc(:))
keeps%nenc(:) = unpack(inserts%nenc(:), lfill_list(:), keeps%nenc(:))

keeps%tpenc1P(:) = unpack(keeps%tpenc1P(:), .not.lfill_list(:), keeps%tpenc1P(:))
keeps%tpenc1P(:) = unpack(inserts%tpenc1P(:), lfill_list(:), keeps%tpenc1P(:))

call whm_fill_pl(keeps, inserts, lfill_list)
class default
write(*,*) 'Error! spill method called for incompatible return type on rmvs_pl'
Expand Down Expand Up @@ -89,12 +84,10 @@ module subroutine rmvs_spill_tp(self, discards, lspill_list)
discards%lperi(:) = pack(keeps%lperi(:), lspill_list(:))
discards%plperP(:) = pack(keeps%plperP(:), lspill_list(:))
discards%plencP(:) = pack(keeps%plencP(:), lspill_list(:))
discards%tpencP(:) = pack(keeps%tpencP(:), lspill_list(:))
if (count(.not.lspill_list(:)) > 0) then
keeps%lperi(:) = pack(keeps%lperi(:), .not. lspill_list(:))
keeps%plperP(:) = pack(keeps%plperP(:), .not. lspill_list(:))
keeps%plencP(:) = pack(keeps%plencP(:), .not. lspill_list(:))
keeps%tpencP(:) = pack(keeps%tpencP(:), .not. lspill_list(:))
end if

call util_spill_tp(keeps, discards, lspill_list)
Expand Down Expand Up @@ -132,9 +125,6 @@ module subroutine rmvs_fill_tp(self, inserts, lfill_list)
keeps%plencP(:) = unpack(keeps%plencP(:), .not.lfill_list(:), keeps%plencP(:))
keeps%plencP(:) = unpack(inserts%plencP(:), lfill_list(:), keeps%plencP(:))

keeps%tpencP(:) = unpack(keeps%tpencP(:), .not.lfill_list(:), keeps%tpencP(:))
keeps%tpencP(:) = unpack(inserts%tpencP(:), lfill_list(:), keeps%tpencP(:))

call util_fill_tp(keeps, inserts, lfill_list)
class default
write(*,*) 'Error! fill method called for incompatible return type on rmvs_tp'
Expand Down
75 changes: 35 additions & 40 deletions src/rmvs/rmvs_step.f90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -70,7 +70,7 @@ module subroutine rmvs_step_out(self, cb, tp, dt, config)
real(DP), intent(in) :: dt !! Step size
class(swiftest_configuration), intent(in) :: config !! Input collection of configuration parameters
! Internals
integer(I4B) :: i, j, k, nenc, itpp
integer(I4B) :: i, j, k
real(DP) :: dto, time

! executable code
Expand Down Expand Up @@ -161,48 +161,46 @@ module subroutine rmvs_step_in_pl(self, cb, tp, config, t, dt)
real(DP), intent(in) :: dt !! Step size
! Internals
logical :: lfirsttp
integer(I4B) :: i, j, nenc
integer(I4B) :: i, j
real(DP) :: dti, time
real(DP), dimension(NDIM, self%nbody) :: xbeg, xend, vbeg

dti = dt / NTPHENC

associate(pl => self, npl => self%nbody, xht => tp%xh, vht => tp%vh, plind => self%plind)
associate(pl => self, npl => self%nbody, xht => tp%xh, vht => tp%vh, plind => self%plind, nenc => pl%nenc)
if (config%loblatecb) call pl%obl_acc_in(cb)
call pl%make_planetocentric(cb, tp, config)
do i = 1, npl
nenc = pl%nenc(i)
if (nenc > 0) then
if (nenc(i) == 0) cycle
! There are inner encounters with this planet...switch to planetocentric coordinates to proceed
time = t
call pl%tpenc(i)%peri_pass(cb, pl, time, dti, .true., 0, nenc, i, config)
! now step the encountering test particles fully through the inner encounter
lfirsttp = .true.
associate(index => pl%tpenc(i)%index, &
xpc => self%tpenc(i)%xh, vpc => self%tpenc(i)%vh, apc => self%tpenc(i)%ah)
pl%tpenc(i)%lfirst = .true.
do index = 1, NTPHENC ! Integrate over the encounter region, using the "substitute" planetocentric systems at each level
xbeg(:,1) = cb%xin(:) - pl%xin(:, i, index - 1)
xend(:,1) = cb%xin(:) - pl%xin(:, i, index)
vbeg(:,1) = cb%vin(:) - pl%vin(:, i, index - 1)
do j = 1, NDIM
xbeg(j,2:npl) = pl%xin(j, plind(i,2:npl), index - 1) - pl%xin(j, i, index - 1)
xend(j,2:npl) = pl%xin(j, plind(i,2:npl), index) - pl%xin(j, i, index)
vbeg(j,2:npl) = pl%vin(j, plind(i,2:npl), index - 1) - pl%vin(j, i, index - 1)
end do
pl%plenc(i)%xh(:,:) = xbeg(:,:)
pl%plenc(i)%vh(:,:) = vbeg(:,:)
call pl%tpenc(i)%set_beg_end(xbeg = xbeg, xend = xend)
call pl%tpenc(i)%step(pl%cbenc(i), pl%plenc(i), config, time, dti)
do j = 1, NDIM
pl%tpenc(i)%xheliocen(j, :) = pl%tpenc(i)%xh(j, :) + pl%xin(j, i, index)
end do
time = config%t + j * dti
call pl%tpenc(i)%peri_pass(cb, pl, time, dti, .false., index, nenc, i, config)
time = t
call pl%tpenc(i)%peri_pass(cb, pl, time, dti, .true., 0, i, config)
! now step the encountering test particles fully through the inner encounter
lfirsttp = .true.
associate(index => pl%tpenc(i)%index, &
xpc => self%tpenc(i)%xh, vpc => self%tpenc(i)%vh, apc => self%tpenc(i)%ah)
pl%tpenc(i)%lfirst = .true.
do index = 1, NTPHENC ! Integrate over the encounter region, using the "substitute" planetocentric systems at each level
xbeg(:,1) = cb%xin(:) - pl%xin(:, i, index - 1)
xend(:,1) = cb%xin(:) - pl%xin(:, i, index)
vbeg(:,1) = cb%vin(:) - pl%vin(:, i, index - 1)
do j = 1, NDIM
xbeg(j,2:npl) = pl%xin(j, plind(i,2:npl), index - 1) - pl%xin(j, i, index - 1)
xend(j,2:npl) = pl%xin(j, plind(i,2:npl), index) - pl%xin(j, i, index)
vbeg(j,2:npl) = pl%vin(j, plind(i,2:npl), index - 1) - pl%vin(j, i, index - 1)
end do
where(pl%tpenc(i)%status(:) == ACTIVE) pl%tpenc(i)%status(:) = INACTIVE
end associate
end if
pl%plenc(i)%xh(:,:) = xbeg(:,:)
pl%plenc(i)%vh(:,:) = vbeg(:,:)
call pl%tpenc(i)%set_beg_end(xbeg = xbeg, xend = xend)
call pl%tpenc(i)%step(pl%cbenc(i), pl%plenc(i), config, time, dti)
do j = 1, NDIM
pl%tpenc(i)%xheliocen(j, :) = pl%tpenc(i)%xh(j, :) + pl%xin(j, i, index)
end do
time = config%t + j * dti
call pl%tpenc(i)%peri_pass(cb, pl, time, dti, .false., index, i, config)
end do
where(pl%tpenc(i)%status(:) == ACTIVE) pl%tpenc(i)%status(:) = INACTIVE
end associate
end do
call pl%end_planetocentric(cb,tp)

Expand All @@ -225,12 +223,8 @@ module subroutine rmvs_step_make_planetocentric(self, cb, tp, config)
class(rmvs_tp), intent(inout) :: tp !! RMVS test particle object
class(swiftest_configuration), intent(in) :: config !! Input collection of configuration parameters
! Internals
integer(I4B) :: i, j, k
type(rmvs_pl) :: cb_as_pl
logical, dimension(:), allocatable :: copyflag

integer(I4B) :: i, j

!allocate(copyflag(self%nbody))
associate(pl => self, npl => self%nbody, nenc => self%nenc, tpenc => self%tpenc, cbenc => self%cbenc, &
plenc => self%plenc, GMpl => self%Gmass)

Expand All @@ -257,6 +251,8 @@ module subroutine rmvs_step_make_planetocentric(self, cb, tp, config)
allocate(tpenc(i)%peri(nenc(i)))
tpenc(i)%status(:) = ACTIVE
tpenc(i)%lperi(:) = .false.
tpenc(i)%plperP(:) = 0
tpenc(i)%name(:) = pack(tp%name(:), pl%encmask(:,i))
!call tp%spill(tpenc(i), pl%encmask(:,i))
! Grab all the encountering test particles and convert them to a planetocentric frame
do j = 1, NDIM
Expand All @@ -277,6 +273,7 @@ module subroutine rmvs_step_make_planetocentric(self, cb, tp, config)
end if
end do
end associate
return
end subroutine rmvs_step_make_planetocentric

module subroutine rmvs_step_end_planetocentric(self, cb, tp)
Expand Down Expand Up @@ -324,9 +321,7 @@ module subroutine rmvs_step_end_planetocentric(self, cb, tp)
end associate
deallocate(self%encmask)


return

end subroutine rmvs_step_end_planetocentric

end submodule s_rmvs_step

0 comments on commit e59535c

Please sign in to comment.