From 917580061c7fa100b5006412c42852aad4d48c91 Mon Sep 17 00:00:00 2001 From: David A Minton Date: Sat, 11 Dec 2021 11:04:16 -0500 Subject: [PATCH] Fixed wrong order of gr_pos_kick operator during Helio steps --- src/helio/helio_step.f90 | 4 ++-- src/symba/symba_gr.f90 | 4 ++-- src/symba/symba_step.f90 | 15 +++++++++------ 3 files changed, 13 insertions(+), 10 deletions(-) diff --git a/src/helio/helio_step.f90 b/src/helio/helio_step.f90 index fb14c9e75..9f68ce37c 100644 --- a/src/helio/helio_step.f90 +++ b/src/helio/helio_step.f90 @@ -56,8 +56,8 @@ module subroutine helio_step_pl(self, system, param, t, dt) call pl%kick(system, param, t, dth, lbeg=.true.) if (param%lgr) call pl%gr_pos_kick(system, param, dth) call pl%drift(system, param, dt) - call pl%kick(system, param, t + dt, dth, lbeg=.false.) if (param%lgr) call pl%gr_pos_kick(system, param, dth) + call pl%kick(system, param, t + dt, dth, lbeg=.false.) call pl%lindrift(cb, dth, lbeg=.false.) call pl%vb2vh(cb) end select @@ -99,8 +99,8 @@ module subroutine helio_step_tp(self, system, param, t, dt) call tp%kick(system, param, t, dth, lbeg=.true.) if (param%lgr) call tp%gr_pos_kick(system, param, dth) call tp%drift(system, param, dt) - call tp%kick(system, param, t + dt, dth, lbeg=.false.) if (param%lgr) call tp%gr_pos_kick(system, param, dth) + call tp%kick(system, param, t + dt, dth, lbeg=.false.) call tp%lindrift(cb, dth, lbeg=.false.) call tp%vb2vh(vbcb = -cb%ptend) end select diff --git a/src/symba/symba_gr.f90 b/src/symba/symba_gr.f90 index 1b0246aa8..5f8b9c5c1 100644 --- a/src/symba/symba_gr.f90 +++ b/src/symba/symba_gr.f90 @@ -23,7 +23,7 @@ module pure subroutine symba_gr_p4_pl(self, system, param, dt) associate(pl => self, npl => self%nbody) select type(system) class is (symba_nbody_system) - do concurrent(i = 1:npl, pl%lmask(i) .and. pl%levelg(i) == system%irec ) + do concurrent(i = 1:npl, pl%lmask(i) .and. (pl%levelg(i) == system%irec) ) call gr_p4_pos_kick(param, pl%xh(:, i), pl%vb(:, i), dt) end do end select @@ -54,7 +54,7 @@ module pure subroutine symba_gr_p4_tp(self, system, param, dt) associate(tp => self, ntp => self%nbody) select type(system) class is (symba_nbody_system) - do concurrent(i = 1:ntp, tp%lmask(i) .and. tp%levelg(i) == system%irec) + do concurrent(i = 1:ntp, tp%lmask(i) .and. (tp%levelg(i) == system%irec)) call gr_p4_pos_kick(param, tp%xh(:, i), tp%vb(:, i), dt) end do end select diff --git a/src/symba/symba_step.f90 b/src/symba/symba_step.f90 index 25616b3d7..983db28f5 100644 --- a/src/symba/symba_step.f90 +++ b/src/symba/symba_step.f90 @@ -83,13 +83,13 @@ module subroutine symba_step_interp_system(self, param, t, dt) call system%recursive_step(param, t, 0) - call pl%kick(system, param, t, dth, lbeg=.false.) if (param%lgr) call pl%gr_pos_kick(system, param, dth) + call pl%kick(system, param, t, dth, lbeg=.false.) call pl%lindrift(cb, dth, lbeg=.false.) call pl%vb2vh(cb) - call tp%kick(system, param, t, dth, lbeg=.false.) if (param%lgr) call tp%gr_pos_kick(system, param, dth) + call tp%kick(system, param, t, dth, lbeg=.false.) call tp%lindrift(cb, dth, lbeg=.false.) call tp%vb2vh(vbcb = -cb%ptend) end select @@ -186,26 +186,29 @@ module recursive subroutine symba_step_recur_system(self, param, t, ireci) call plplenc_list%kick(system, dth, irecp, -1) call pltpenc_list%kick(system, dth, irecp, -1) end if + if (param%lgr) then call pl%gr_pos_kick(system, param, dth) call tp%gr_pos_kick(system, param, dth) end if + call pl%drift(system, param, dtl) call tp%drift(system, param, dtl) if (lencounter) call system%recursive_step(param, t+dth,irecp) system%irec = ireci + if (param%lgr) then + call pl%gr_pos_kick(system, param, dth) + call tp%gr_pos_kick(system, param, dth) + end if + call plplenc_list%kick(system, dth, irecp, 1) call pltpenc_list%kick(system, dth, irecp, 1) if (ireci /= 0) then call plplenc_list%kick(system, dth, irecp, -1) call pltpenc_list%kick(system, dth, irecp, -1) end if - if (param%lgr) then - call pl%gr_pos_kick(system, param, dth) - call tp%gr_pos_kick(system, param, dth) - end if if (param%lclose) then lplpl_collision = plplenc_list%collision_check(system, param, t+dtl, dtl, ireci)