diff --git a/src/modules/swiftest_classes.f90 b/src/modules/swiftest_classes.f90 index f255e1a38..aed6f23b0 100644 --- a/src/modules/swiftest_classes.f90 +++ b/src/modules/swiftest_classes.f90 @@ -260,8 +260,9 @@ module swiftest_classes procedure :: append => util_append_tp !! Appends elements from one structure to another procedure :: h2b => util_coord_h2b_tp !! Convert test particles from heliocentric to barycentric coordinates (position and velocity) procedure :: b2h => util_coord_b2h_tp !! Convert test particles from barycentric to heliocentric coordinates (position and velocity) - procedure :: vh2vb => util_coord_vh2vb_tp !! Convert test particles from heliocentric to barycentric coordinates (velocity only) procedure :: vb2vh => util_coord_vb2vh_tp !! Convert test particles from barycentric to heliocentric coordinates (velocity only) + procedure :: vh2vb => util_coord_vh2vb_tp !! Convert test particles from heliocentric to barycentric coordinates (velocity only) + procedure :: xh2xb => util_coord_xh2xb_tp !! Convert test particles from heliocentric to barycentric coordinates (position only) procedure :: fill => util_fill_tp !! "Fills" bodies from one object into another depending on the results of a mask (uses the UNPACK intrinsic) procedure :: get_peri => util_peri_tp !! Determine system pericenter passages for test particles procedure :: resize => util_resize_tp !! Checks the current size of a Swiftest body against the requested size and resizes it if it is too small. diff --git a/src/modules/swiftest_globals.f90 b/src/modules/swiftest_globals.f90 index cbe626e43..454f454f5 100644 --- a/src/modules/swiftest_globals.f90 +++ b/src/modules/swiftest_globals.f90 @@ -72,23 +72,24 @@ module swiftest_globals integer(I4B), parameter :: HYPERBOLA = 1 !! Symbolic names for orbit types - hyperbola !> Symbolic names for body/particle status codes: - integer(I4B), parameter :: ACTIVE = 0 - integer(I4B), parameter :: INACTIVE = 1 - integer(I4B), parameter :: DISCARDED_RMAX = -1 - integer(I4B), parameter :: DISCARDED_RMIN = -2 - integer(I4B), parameter :: DISCARDED_RMAXU = -3 - integer(I4B), parameter :: DISCARDED_PERI = -4 - integer(I4B), parameter :: DISCARDED_PLR = -5 - integer(I4B), parameter :: DISCARDED_PLQ = -6 - integer(I4B), parameter :: DISCARDED_DRIFTERR = -7 - integer(I4B), parameter :: MERGED = -8 - integer(I4B), parameter :: DISRUPTION = -9 - integer(I4B), parameter :: SUPERCATASTROPHIC = -10 - integer(I4B), parameter :: GRAZE_AND_MERGE = -11 - integer(I4B), parameter :: HIT_AND_RUN = -12 - integer(I4B), parameter :: COLLISION = -13 - integer(I4B), parameter :: NEW_PARTICLE = -14 - integer(I4B), parameter :: OLD_PARTICLE = -15 + integer(I4B), parameter :: ACTIVE = 0 + integer(I4B), parameter :: INACTIVE = 1 + integer(I4B), parameter :: DISCARDED_RMAX = -1 + integer(I4B), parameter :: DISCARDED_RMIN = -2 + integer(I4B), parameter :: DISCARDED_RMAXU = -3 + integer(I4B), parameter :: DISCARDED_PERI = -4 + integer(I4B), parameter :: DISCARDED_PLR = -5 + integer(I4B), parameter :: DISCARDED_PLQ = -6 + integer(I4B), parameter :: DISCARDED_DRIFTERR = -7 + integer(I4B), parameter :: MERGED = -8 + integer(I4B), parameter :: DISRUPTION = -9 + integer(I4B), parameter :: SUPERCATASTROPHIC = -10 + integer(I4B), parameter :: GRAZE_AND_MERGE = -11 + integer(I4B), parameter :: HIT_AND_RUN_DISRUPT = -12 + integer(I4B), parameter :: HIT_AND_RUN_PURE = -13 + integer(I4B), parameter :: COLLISION = -14 + integer(I4B), parameter :: NEW_PARTICLE = -15 + integer(I4B), parameter :: OLD_PARTICLE = -16 !>Symbolic names for collisional outcomes from collresolve_resolve: integer(I4B), parameter :: COLLRESOLVE_REGIME_MERGE = 1 diff --git a/src/symba/symba_collision.f90 b/src/symba/symba_collision.f90 index 2f93881b1..1a13bdfd0 100644 --- a/src/symba/symba_collision.f90 +++ b/src/symba/symba_collision.f90 @@ -186,7 +186,7 @@ module function symba_collision_casehitandrun(system, param, family, x, v, mass, end if end if if (lpure) then ! Reset these bodies back to being active so that nothing further is done to them - status = ACTIVE + status = HIT_AND_RUN_PURE select type(pl => system%pl) class is (symba_pl) pl%status(family(:)) = status @@ -194,7 +194,7 @@ module function symba_collision_casehitandrun(system, param, family, x, v, mass, pl%lcollision(family(:)) = .false. end select else - status = HIT_AND_RUN + status = HIT_AND_RUN_DISRUPT call symba_collision_mergeaddsub(system, param, family, id_frag, Ip_frag, m_frag, rad_frag, xb_frag, vb_frag, rot_frag, status) end if @@ -871,7 +871,7 @@ subroutine symba_collision_mergeaddsub(system, param, family, id_frag, Ip_frag, plnew%info(i)%origin_xh(:) = plnew%xh(:,i) plnew%info(i)%origin_vh(:) = plnew%vh(:,i) end do - case(HIT_AND_RUN) + case(HIT_AND_RUN_DISRUPT) plnew%info(1) = pl%info(ibiggest) plnew%status(1) = OLD_PARTICLE plnew%status(2:nfrag) = NEW_PARTICLE @@ -938,7 +938,7 @@ module subroutine symba_collision_resolve_fragmentations(self, system, param) real(DP), dimension(NDIM,2) :: x, v, L_spin, Ip !! Output values that represent a 2-body equivalent of a possibly 2+ body collision real(DP), dimension(2) :: mass, radius !! Output values that represent a 2-body equivalent of a possibly 2+ body collision logical :: lgoodcollision - integer(I4B) :: i, status, jtarg, jproj, regime + integer(I4B) :: i, jtarg, jproj, regime real(DP), dimension(2) :: radius_si, mass_si, density_si real(DP) :: mtiny_si, Mcb_si real(DP), dimension(NDIM) :: x1_si, v1_si, x2_si, v2_si @@ -994,13 +994,13 @@ module subroutine symba_collision_resolve_fragmentations(self, system, param) select case (regime) case (COLLRESOLVE_REGIME_DISRUPTION) - status = symba_collision_casedisruption(system, param, family, x, v, mass, radius, L_spin, Ip, mass_res, Qloss) + plpl_collisions%status = symba_collision_casedisruption(system, param, family, x, v, mass, radius, L_spin, Ip, mass_res, Qloss) case (COLLRESOLVE_REGIME_SUPERCATASTROPHIC) - status = symba_collision_casesupercatastrophic(system, param, family, x, v, mass, radius, L_spin, Ip, mass_res, Qloss) + plpl_collisions%status = symba_collision_casesupercatastrophic(system, param, family, x, v, mass, radius, L_spin, Ip, mass_res, Qloss) case (COLLRESOLVE_REGIME_HIT_AND_RUN) - status = symba_collision_casehitandrun(system, param, family, x, v, mass, radius, L_spin, Ip, mass_res, Qloss) + plpl_collisions%status = symba_collision_casehitandrun(system, param, family, x, v, mass, radius, L_spin, Ip, mass_res, Qloss) case (COLLRESOLVE_REGIME_MERGE, COLLRESOLVE_REGIME_GRAZE_AND_MERGE) - status = symba_collision_casemerge(system, param, family, x, v, mass, radius, L_spin, Ip) + plpl_collisions%status = symba_collision_casemerge(system, param, family, x, v, mass, radius, L_spin, Ip) case default write(*,*) "Error in symba_collision, unrecognized collision regime" call util_exit(FAILURE) @@ -1140,12 +1140,11 @@ 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) call system%tp%discard(system, param) return end subroutine symba_collision_resolve_pltpenc - - end submodule s_symba_collision \ No newline at end of file