From 53ca9640bc833f5d815ffe99e6eb7fe869a60747 Mon Sep 17 00:00:00 2001 From: David Minton Date: Fri, 2 Apr 2021 23:21:19 -0400 Subject: [PATCH] Improved encounter check memory requirements --- .../mars_ejecta/config.swiftest.in | 2 +- .../mars_ejecta/param.swifter.in | 2 +- .../mars_ejecta/profswifter.sh | 2 ++ .../rmvs_swifter_comparison/mars_ejecta/tp.in | 2 +- src/rmvs/rmvs_encounter_check.f90 | 17 ++++++++--------- src/rmvs/rmvs_step.f90 | 19 +++++++++++++------ 6 files changed, 26 insertions(+), 18 deletions(-) create mode 100755 examples/rmvs_swifter_comparison/mars_ejecta/profswifter.sh diff --git a/examples/rmvs_swifter_comparison/mars_ejecta/config.swiftest.in b/examples/rmvs_swifter_comparison/mars_ejecta/config.swiftest.in index 685ee4595..fd124e6a5 100644 --- a/examples/rmvs_swifter_comparison/mars_ejecta/config.swiftest.in +++ b/examples/rmvs_swifter_comparison/mars_ejecta/config.swiftest.in @@ -4,7 +4,7 @@ !NPLMAX -1 ! not used !NTPMAX -1 ! not used T0 0.0e0 -TSTOP 365.25e5 ! simulation length in days +TSTOP 365.25e3 ! simulation length in days DT 2e0 ! stepsize in days ISTEP_OUT 73500 ! output cadence ISTEP_DUMP 73500 ! system dump cadence diff --git a/examples/rmvs_swifter_comparison/mars_ejecta/param.swifter.in b/examples/rmvs_swifter_comparison/mars_ejecta/param.swifter.in index 5327ebcab..67e4ed92a 100644 --- a/examples/rmvs_swifter_comparison/mars_ejecta/param.swifter.in +++ b/examples/rmvs_swifter_comparison/mars_ejecta/param.swifter.in @@ -4,7 +4,7 @@ !NPLMAX -1 ! not used !NTPMAX -1 ! not used T0 0.0e0 -TSTOP 365.25e5 ! simulation length in days +TSTOP 365.25e3 ! simulation length in days DT 2e0 ! stepsize in days ISTEP_OUT 73500 ! output cadence ISTEP_DUMP 73500 ! system dump cadence diff --git a/examples/rmvs_swifter_comparison/mars_ejecta/profswifter.sh b/examples/rmvs_swifter_comparison/mars_ejecta/profswifter.sh new file mode 100755 index 000000000..a03493f54 --- /dev/null +++ b/examples/rmvs_swifter_comparison/mars_ejecta/profswifter.sh @@ -0,0 +1,2 @@ +#!/bin/bash +gprof ./swifter_rmvs | /home/daminton/git/gprof2dot/gprof2dot.py | dot -Tpng -o swifter_profile.png diff --git a/examples/rmvs_swifter_comparison/mars_ejecta/tp.in b/examples/rmvs_swifter_comparison/mars_ejecta/tp.in index f2fb549d8..d5eb0e3f8 100644 --- a/examples/rmvs_swifter_comparison/mars_ejecta/tp.in +++ b/examples/rmvs_swifter_comparison/mars_ejecta/tp.in @@ -1,4 +1,4 @@ - 10 + 200 6000200 0.321794727714005 -1.39371227734394 -3.555372224179648E-002 1.263074812292130E-002 3.554726326213714E-003 5.297001107625803E-004 diff --git a/src/rmvs/rmvs_encounter_check.f90 b/src/rmvs/rmvs_encounter_check.f90 index 63ead1338..030d181c9 100644 --- a/src/rmvs/rmvs_encounter_check.f90 +++ b/src/rmvs/rmvs_encounter_check.f90 @@ -21,26 +21,25 @@ module function rmvs_encounter_check_tp(self, cb, pl, dt, rts) result(lencounter real(DP) :: r2, v2, vdotr real(DP), dimension(NDIM) :: xr, vr real(DP), dimension(pl%nbody) :: r2crit + logical :: lflag associate(tp => self, ntp => self%nbody, npl => pl%nbody, rhill => pl%rhill, & xht => self%xh, vht => self%vh, xbeg => self%xbeg, vbeg => self%vbeg) - if (.not.allocated(pl%encmask)) allocate(pl%encmask(ntp, npl)) - pl%encmask(:,:) = .false. r2crit(:) = (rts * rhill(:))**2 - do i = 1, ntp - if (tp%status(i) /= ACTIVE) cycle - tp%plencP(i) = 0 - do j = 1, npl + tp%plencP(:) = 0 + do j = 1, npl + do i = 1, ntp + if ((tp%status(i) /= ACTIVE).or.(tp%plencP(i) /=0)) cycle xr(:) = xht(:, i) - xbeg(:, j) vr(:) = vht(:, i) - vbeg(:, j) r2 = dot_product(xr(:), xr(:)) v2 = dot_product(vr(:), vr(:)) vdotr = dot_product(vr(:), xr(:)) - pl%encmask(i,j) = rmvs_chk_ind(r2, v2, vdotr, dt, r2crit(j)) - if (pl%encmask(i,j)) tp%plencP(i) = j + lflag = rmvs_chk_ind(r2, v2, vdotr, dt, r2crit(j)) + if (lflag) tp%plencP(i) = j end do + pl%nenc(j) = count(tp%plencP(:) == j) end do - pl%nenc(:) = count(pl%encmask(:,:), dim=1) lencounter = any(pl%nenc(:) > 0) end associate return diff --git a/src/rmvs/rmvs_step.f90 b/src/rmvs/rmvs_step.f90 index 97b2303be..b665ece45 100644 --- a/src/rmvs/rmvs_step.f90 +++ b/src/rmvs/rmvs_step.f90 @@ -224,6 +224,7 @@ module subroutine rmvs_step_make_planetocentric(self, cb, tp, config) class(swiftest_configuration), intent(in) :: config !! Input collection of configuration parameters ! Internals integer(I4B) :: i, j + logical, dimension(:), allocatable :: encmask associate(pl => self, npl => self%nbody, nenc => self%nenc, tpenc => self%tpenc, cbenc => self%cbenc, & plenc => self%plenc, GMpl => self%Gmass) @@ -231,6 +232,9 @@ module subroutine rmvs_step_make_planetocentric(self, cb, tp, config) do i = 1, npl if (nenc(i) == 0) cycle ! There are inner encounters with this planet + if (allocated(encmask)) deallocate(encmask) + allocate(encmask(nenc(i))) + encmask(:) = tp%plencP(:) == i ! Save the index value of the planet corresponding to this encounter tpenc(i)%ipleP = i @@ -252,13 +256,13 @@ module subroutine rmvs_step_make_planetocentric(self, cb, tp, config) tpenc(i)%status(:) = ACTIVE tpenc(i)%lperi(:) = .false. tpenc(i)%plperP(:) = 0 - tpenc(i)%name(:) = pack(tp%name(:), pl%encmask(:,i)) + tpenc(i)%name(:) = pack(tp%name(:), encmask(:)) !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 - tpenc(i)%xheliocen(j, :) = pack(tp%xh(j,:), pl%encmask(:,i)) + tpenc(i)%xheliocen(j, :) = pack(tp%xh(j,:), encmask(:)) tpenc(i)%xh(j, :) = tpenc(i)%xheliocen(j, :) - pl%xin(j, i, 0) - tpenc(i)%vh(j, :) = pack(tp%vh(j,:), pl%encmask(:,i)) - pl%vin(j, i, 0) + tpenc(i)%vh(j, :) = pack(tp%vh(j,:), encmask(:)) - pl%vin(j, i, 0) end do ! Make sure that the test particles get the planetocentric value of mu @@ -289,15 +293,19 @@ module subroutine rmvs_step_end_planetocentric(self, cb, tp) ! Internals integer(I4B) :: i, j integer(I4B), dimension(:), allocatable :: tpind + logical, dimension(:), allocatable :: encmask associate(pl => self, nenc => self%nenc, npl => self%nbody, ntp => tp%nbody, name => tp%name, & - tpenc => self%tpenc, plenc => self%plenc, cbenc => self%cbenc, encmask => self%encmask) + tpenc => self%tpenc, plenc => self%plenc, cbenc => self%cbenc) do i = 1, npl if (nenc(i) == 0) cycle allocate(tpind(nenc(i))) ! Index array of encountering test particles - tpind(:) = pack([(j,j=1,ntp)], encmask(1:ntp, i)) + if (allocated(encmask)) deallocate(encmask) + allocate(encmask(nenc(i))) + encmask(:) = tp%plencP(:) == i + tpind(:) = pack([(j,j=1,ntp)], encmask(:)) ! Copy the results of the integration back over and shift back to heliocentric reference tp%status(tpind(1:nenc(i))) = tpenc(i)%status(1:nenc(i)) @@ -319,7 +327,6 @@ module subroutine rmvs_step_end_planetocentric(self, cb, tp) deallocate(tpind) end do end associate - deallocate(self%encmask) return end subroutine rmvs_step_end_planetocentric