Made sure to sync up the coordinate systems before doing a tp discard…

…. Otherwise the h2b call in discard.f90 uses the pre-encounter planet heliocentric position and velocity vectors to compute barycentric, instead of using the current barycentric vectors when doing a SyMBA recursion step.

src/orbel/orbel.f90

@@ -883,6 +883,8 @@ module subroutine orbel_xv2el_vec(self, cb)
                           self%a(i), self%e(i), self%inc(i),  &
                           self%capom(i), self%omega(i), self%capm(i))
       end do
+
+      return
    end subroutine orbel_xv2el_vec 
 
 

src/symba/symba_collision.f90

@@ -370,7 +370,7 @@ module function symba_collision_check_encounter(self, system, param, t, dt, irec
          end select
       end select
 
-      ! Extract the pl-pl or pl-tpencounter list and return the plplcollision_list
+      ! Extract the pl-pl or pl-tp encounter list and return the pl-pl or pl-tp collision_list
       if (lany_collision) then
          select type(self)
          class is (symba_plplenc)
@@ -1074,8 +1074,14 @@ module subroutine symba_collision_resolve_pltpenc(self, system, param, t, dt, ir
       real(DP),                   intent(in)    :: t      !! Current simulation tim
       real(DP),                   intent(in)    :: dt     !! Current simulation step size
       integer(I4B),               intent(in)    :: irec   !! Current recursion level
-
-      call system%tp%xh2xb(system%cb)
+
+      ! Make sure coordinate systems are all synced up due to being inside the recursion at this point
+      call system%pl%vb2vh(system%cb)
+      call system%tp%vb2vh(system%cb%vb)
+      call system%pl%b2h(system%cb)
+      call system%tp%b2h(system%cb)
+
+      ! Discard the collider
       call system%tp%discard(system, param)
 
       return