diff --git a/src/collision/collision_io.f90 b/src/collision/collision_io.f90 index 2c5abd7b9..864fa93f9 100644 --- a/src/collision/collision_io.f90 +++ b/src/collision/collision_io.f90 @@ -191,7 +191,7 @@ module subroutine collision_io_initialize_output(self, param) 667 continue write(*,*) "Error creating fragmentation output file. " // trim(adjustl(errmsg)) - call util_exit(FAILURE) + call swiftest_util_exit(FAILURE) end subroutine collision_io_initialize_output diff --git a/src/collision/collision_resolve.f90 b/src/collision/collision_resolve.f90 index 4ac8bdf16..f08a056ac 100644 --- a/src/collision/collision_resolve.f90 +++ b/src/collision/collision_resolve.f90 @@ -635,7 +635,7 @@ subroutine collision_resolve_list(plpl_collision , system, param, t) plpl_collision%status(i) = collision_resolve_merge(system, param, t) case default write(*,*) "Error in swiftest_collision, unrecognized collision regime" - call util_exit(FAILURE) + call swiftest_util_exit(FAILURE) end select call collision_history%take_snapshot(param,system, t, "after") call impactors%reset() diff --git a/src/collision/collision_util.f90 b/src/collision/collision_util.f90 index 0690adf60..02eaafdd3 100644 --- a/src/collision/collision_util.f90 +++ b/src/collision/collision_util.f90 @@ -291,7 +291,7 @@ module subroutine collision_util_get_energy_momentum(self, system, param, lbefo if (.not.allocated(tmpsys)) then write(*,*) "Error in collision_util_get_energy_momentum. " // & " This must be called with lbefore=.true. at least once before calling it with lbefore=.false." - call util_exit(FAILURE) + call swiftest_util_exit(FAILURE) end if select type(tmpsys) class is (swiftest_nbody_system) @@ -345,7 +345,7 @@ module subroutine collision_util_index_map(self) call self%get_index_values(idvals, tvals) ! Consolidate ids to only unique values - call util_unique(idvals,self%idvals,self%idmap) + call swiftest_util_unique(idvals,self%idvals,self%idmap) self%nid = size(self%idvals) ! Don't consolidate time values (multiple collisions can happen in a single time step) @@ -502,6 +502,7 @@ module subroutine collision_util_set_coordinate_system(self) return end subroutine collision_util_set_coordinate_system + subroutine collision_util_save_snapshot(collision_history, snapshot) !! author: David A. Minton !! diff --git a/src/encounter/encounter_check.f90 b/src/encounter/encounter_check.f90 index f3a5a3f8d..9aad9cfdb 100644 --- a/src/encounter/encounter_check.f90 +++ b/src/encounter/encounter_check.f90 @@ -168,10 +168,10 @@ module subroutine encounter_check_all_plplm(param, nplm, nplt, xplm, vplm, xplt, call move_alloc(ltmp, lvdotr) nenc = nenc + plmplt_nenc - call util_sort(index1, ind) - call util_sort_rearrange(index1, ind, nenc) - call util_sort_rearrange(index2, ind, nenc) - call util_sort_rearrange(lvdotr, ind, nenc) + call swiftest_util_sort(index1, ind) + call swiftest_util_sort_rearrange(index1, ind, nenc) + call swiftest_util_sort_rearrange(index2, ind, nenc) + call swiftest_util_sort_rearrange(lvdotr, ind, nenc) end if @@ -565,7 +565,7 @@ subroutine encounter_check_all_triangular_plpl(npl, x, v, renc, dt, nenc, index1 integer(I4B), dimension(:), allocatable, save :: ind_arr type(collision_list_plpl), dimension(npl) :: lenc - call util_index_array(ind_arr, npl) + call swiftest_util_index_array(ind_arr, npl) !$omp parallel do default(private) schedule(static)& !$omp shared(x, v, renc, lenc, ind_arr) & @@ -612,7 +612,7 @@ subroutine encounter_check_all_triangular_plplm(nplm, nplt, xplm, vplm, xplt, vp integer(I4B), dimension(:), allocatable, save :: ind_arr type(collision_list_plpl), dimension(nplm) :: lenc - call util_index_array(ind_arr, nplt) + call swiftest_util_index_array(ind_arr, nplt) !$omp parallel do default(private) schedule(dynamic)& !$omp shared(xplm, vplm, xplt, vplt, rencm, renct, lenc, ind_arr) & @@ -659,7 +659,7 @@ subroutine encounter_check_all_triangular_pltp(npl, ntp, xpl, vpl, xtp, vtp, ren type(collision_list_pltp), dimension(npl) :: lenc real(DP), dimension(ntp) :: renct - call util_index_array(ind_arr, ntp) + call swiftest_util_index_array(ind_arr, ntp) renct(:) = 0.0_DP !$omp parallel do default(private) schedule(dynamic)& @@ -807,10 +807,10 @@ subroutine encounter_check_remove_duplicates(n, nenc, index1, index2, lvdotr) return end if - call util_sort(index1, ind) - call util_sort_rearrange(index1, ind, nenc) - call util_sort_rearrange(index2, ind, nenc) - call util_sort_rearrange(lvdotr, ind, nenc) + call swiftest_util_sort(index1, ind) + call swiftest_util_sort_rearrange(index1, ind, nenc) + call swiftest_util_sort_rearrange(index2, ind, nenc) + call swiftest_util_sort_rearrange(lvdotr, ind, nenc) ! Get the bounds on the bodies in the first index ibeg(:) = n @@ -836,7 +836,7 @@ subroutine encounter_check_remove_duplicates(n, nenc, index1, index2, lvdotr) khi = iend(i) nenci = khi - klo + 1_I8B if (allocated(ind)) deallocate(ind) - call util_sort(index2(klo:khi), ind) + call swiftest_util_sort(index2(klo:khi), ind) index2(klo:khi) = itmp(klo - 1_I8B + ind(:)) do j = klo + 1_I8B, khi if (index2(j) == index2(j - 1_I8B)) lencounter(j) = .false. @@ -876,7 +876,7 @@ pure module subroutine encounter_check_sort_aabb_1D(self, n, extent_arr) ! Internals integer(I8B) :: i, k - call util_sort(extent_arr, self%ind) + call swiftest_util_sort(extent_arr, self%ind) do concurrent(k = 1_I8B:2_I8B * n) i = self%ind(k) @@ -923,7 +923,7 @@ module subroutine encounter_check_sweep_aabb_double_list(self, n1, n2, r1, v1, r real(DP), dimension(2*(n1+n2)) :: xind, yind, zind, vxind, vyind, vzind, rencind ntot = n1 + n2 - call util_index_array(ind_arr, ntot) + call swiftest_util_index_array(ind_arr, ntot) do concurrent(dim = 1:SWEEPDIM) loverlap_by_dimension(dim,:) = (self%aabb(dim)%ibeg(:) + 1_I8B) < (self%aabb(dim)%iend(:) - 1_I8B) @@ -1035,7 +1035,7 @@ module subroutine encounter_check_sweep_aabb_single_list(self, n, x, v, renc, dt integer(I4B), dimension(:), allocatable, save :: ind_arr integer(I8B) :: ibeg, iend - call util_index_array(ind_arr, n) + call swiftest_util_index_array(ind_arr, n) dim = 1 ! Sweep the intervals for each of the massive bodies along one dimension diff --git a/src/encounter/encounter_io.f90 b/src/encounter/encounter_io.f90 index b92b9a2fd..df82d567e 100644 --- a/src/encounter/encounter_io.f90 +++ b/src/encounter/encounter_io.f90 @@ -146,7 +146,7 @@ module subroutine encounter_io_initialize_output(self, param) 667 continue write(*,*) "Error creating encounter output file. " // trim(adjustl(errmsg)) - call util_exit(FAILURE) + call swiftest_util_exit(FAILURE) end subroutine encounter_io_initialize_output diff --git a/src/encounter/encounter_setup.f90 b/src/encounter/encounter_setup.f90 index aff1bd626..3cf165129 100644 --- a/src/encounter/encounter_setup.f90 +++ b/src/encounter/encounter_setup.f90 @@ -84,6 +84,7 @@ module subroutine encounter_setup_list(self, n) allocate(self%r2(NDIM,n)) allocate(self%v1(NDIM,n)) allocate(self%v2(NDIM,n)) + allocate(self%level(n)) self%tcollision(:) = 0.0_DP self%lvdotr(:) = .false. @@ -97,6 +98,7 @@ module subroutine encounter_setup_list(self, n) self%r2(:,:) = 0.0_DP self%v1(:,:) = 0.0_DP self%v2(:,:) = 0.0_DP + self%level(:) = 0 return end subroutine encounter_setup_list diff --git a/src/encounter/encounter_util.f90 b/src/encounter/encounter_util.f90 index 8dc63201e..b124c5fbf 100644 --- a/src/encounter/encounter_util.f90 +++ b/src/encounter/encounter_util.f90 @@ -26,17 +26,19 @@ module subroutine encounter_util_append_list(self, source, lsource_mask) nold = self%nenc nsrc = source%nenc - call util_append(self%lvdotr, source%lvdotr, nold, nsrc, lsource_mask) - call util_append(self%lclosest, source%lclosest, nold, nsrc, lsource_mask) - call util_append(self%status, source%status, nold, nsrc, lsource_mask) - call util_append(self%index1, source%index1, nold, nsrc, lsource_mask) - call util_append(self%index2, source%index2, nold, nsrc, lsource_mask) - call util_append(self%id1, source%id1, nold, nsrc, lsource_mask) - call util_append(self%id2, source%id2, nold, nsrc, lsource_mask) - call util_append(self%r1, source%r1, nold, nsrc, lsource_mask) - call util_append(self%r2, source%r2, nold, nsrc, lsource_mask) - call util_append(self%v1, source%v1, nold, nsrc, lsource_mask) - call util_append(self%v2, source%v2, nold, nsrc, lsource_mask) + call swiftest_util_append(self%tcollision, source%tcollision, nold, nsrc, lsource_mask) + call swiftest_util_append(self%lclosest, source%lclosest, nold, nsrc, lsource_mask) + call swiftest_util_append(self%lvdotr, source%lvdotr, nold, nsrc, lsource_mask) + call swiftest_util_append(self%status, source%status, nold, nsrc, lsource_mask) + call swiftest_util_append(self%index1, source%index1, nold, nsrc, lsource_mask) + call swiftest_util_append(self%index2, source%index2, nold, nsrc, lsource_mask) + call swiftest_util_append(self%id1, source%id1, nold, nsrc, lsource_mask) + call swiftest_util_append(self%id2, source%id2, nold, nsrc, lsource_mask) + call swiftest_util_append(self%r1, source%r1, nold, nsrc, lsource_mask) + call swiftest_util_append(self%r2, source%r2, nold, nsrc, lsource_mask) + call swiftest_util_append(self%v1, source%v1, nold, nsrc, lsource_mask) + call swiftest_util_append(self%v2, source%v2, nold, nsrc, lsource_mask) + call swiftest_util_append(self%level, source%level, nold, nsrc, lsource_mask) self%nenc = nold + count(lsource_mask(1:nsrc)) return @@ -55,8 +57,10 @@ module subroutine encounter_util_copy_list(self, source) associate(n => source%nenc) self%nenc = n self%t = source%t - self%lvdotr(1:n) = source%lvdotr(1:n) + self%lcollision = source%lcollision + self%tcollision(1:n) = source%tcollision(1:n) self%lclosest(1:n) = source%lclosest(1:n) + self%lvdotr(1:n) = source%lvdotr(1:n) self%status(1:n) = source%status(1:n) self%index1(1:n) = source%index1(1:n) self%index2(1:n) = source%index2(1:n) @@ -66,6 +70,7 @@ module subroutine encounter_util_copy_list(self, source) self%r2(:,1:n) = source%r2(:,1:n) self%v1(:,1:n) = source%v1(:,1:n) self%v2(:,1:n) = source%v2(:,1:n) + self%level(1:n) = source%level(1:n) end associate return @@ -97,8 +102,8 @@ module subroutine encounter_util_dealloc_list(self) class(encounter_list), intent(inout) :: self if (allocated(self%tcollision)) deallocate(self%tcollision) - if (allocated(self%lvdotr)) deallocate(self%lvdotr) if (allocated(self%lclosest)) deallocate(self%lclosest) + if (allocated(self%lvdotr)) deallocate(self%lvdotr) if (allocated(self%status)) deallocate(self%status) if (allocated(self%index1)) deallocate(self%index1) if (allocated(self%index2)) deallocate(self%index2) @@ -108,6 +113,7 @@ module subroutine encounter_util_dealloc_list(self) if (allocated(self%r2)) deallocate(self%r2) if (allocated(self%v1)) deallocate(self%v1) if (allocated(self%v2)) deallocate(self%v2) + if (allocated(self%level)) deallocate(self%level) return end subroutine encounter_util_dealloc_list @@ -276,11 +282,11 @@ module subroutine encounter_util_index_map(self) call encounter_util_get_vals_storage(self, idvals, tvals) ! Consolidate ids to only unique values - call util_unique(idvals,self%idvals,self%idmap) + call swiftest_util_unique(idvals,self%idvals,self%idmap) self%nid = size(self%idvals) ! Consolidate time values to only unique values - call util_unique(tvals,self%tvals,self%tmap) + call swiftest_util_unique(tvals,self%tvals,self%tmap) self%nt = size(self%tvals) return @@ -339,17 +345,19 @@ module subroutine encounter_util_spill_list(self, discards, lspill_list, ldestru integer(I8B) :: nenc_old associate(keeps => self) - call util_spill(keeps%lvdotr, discards%lvdotr, lspill_list, ldestructive) - call util_spill(keeps%lclosest, discards%lclosest, lspill_list, ldestructive) - call util_spill(keeps%status, discards%status, lspill_list, ldestructive) - call util_spill(keeps%index1, discards%index1, lspill_list, ldestructive) - call util_spill(keeps%index2, discards%index2, lspill_list, ldestructive) - call util_spill(keeps%id1, discards%id1, lspill_list, ldestructive) - call util_spill(keeps%id2, discards%id2, lspill_list, ldestructive) - call util_spill(keeps%r1, discards%r1, lspill_list, ldestructive) - call util_spill(keeps%r2, discards%r2, lspill_list, ldestructive) - call util_spill(keeps%v1, discards%v1, lspill_list, ldestructive) - call util_spill(keeps%v2, discards%v2, lspill_list, ldestructive) + call swiftest_util_spill(keeps%tcollision, discards%tcollision, lspill_list, ldestructive) + call swiftest_util_spill(keeps%lvdotr, discards%lvdotr, lspill_list, ldestructive) + call swiftest_util_spill(keeps%lclosest, discards%lclosest, lspill_list, ldestructive) + call swiftest_util_spill(keeps%status, discards%status, lspill_list, ldestructive) + call swiftest_util_spill(keeps%index1, discards%index1, lspill_list, ldestructive) + call swiftest_util_spill(keeps%index2, discards%index2, lspill_list, ldestructive) + call swiftest_util_spill(keeps%id1, discards%id1, lspill_list, ldestructive) + call swiftest_util_spill(keeps%id2, discards%id2, lspill_list, ldestructive) + call swiftest_util_spill(keeps%r1, discards%r1, lspill_list, ldestructive) + call swiftest_util_spill(keeps%r2, discards%r2, lspill_list, ldestructive) + call swiftest_util_spill(keeps%v1, discards%v1, lspill_list, ldestructive) + call swiftest_util_spill(keeps%v2, discards%v2, lspill_list, ldestructive) + call swiftest_util_spill(keeps%level, discards%level, lspill_list, ldestructive) nenc_old = keeps%nenc diff --git a/src/fraggle/fraggle_generate.f90 b/src/fraggle/fraggle_generate.f90 index 5f9eccee8..59397a808 100644 --- a/src/fraggle/fraggle_generate.f90 +++ b/src/fraggle/fraggle_generate.f90 @@ -401,7 +401,7 @@ subroutine fraggle_generate_tan_vel(collision_system, lfailure) tol = TOL_INIT do while(tol < TOL_MIN) - call util_minimize_bfgs(objective_function, nfrag-6, v_t_initial(7:nfrag), tol, MAXLOOP, lfailure, v_t_output) + call swiftest_util_minimize_bfgs(objective_function, nfrag-6, v_t_initial(7:nfrag), tol, MAXLOOP, lfailure, v_t_output) fragments%v_t_mag(7:nfrag) = v_t_output(:) ! Now that the KE-minimized values of the i>6 fragments are found, calculate the momentum-conserving solution for tangential velociteis v_t_initial(7:nfrag) = fragments%v_t_mag(7:nfrag) @@ -495,7 +495,7 @@ function solve_fragment_tan_vel(lfailure, v_t_mag_input) result(v_t_mag_output) b(1:3) = -L_lin_others(:) b(4:6) = fragments%L_budget(:) - fragments%Lspin(:) - L_orb_others(:) allocate(v_t_mag_output(nfrag)) - v_t_mag_output(1:6) = util_solve_linear_system(A, b, 6, lfailure) + v_t_mag_output(1:6) = swiftest_util_solve_linear_system(A, b, 6, lfailure) if (present(v_t_mag_input)) v_t_mag_output(7:nfrag) = v_t_mag_input(:) end associate end select @@ -559,7 +559,7 @@ subroutine fraggle_generate_rad_vel(collision_system, lfailure) integer(I4B), parameter :: MAXLOOP = 100 real(DP) :: ke_radial, tol integer(I4B) :: i - real(DP), dimension(:), allocatable :: v_r_initial + real(DP), dimension(:), allocatable :: v_r_initial, v_r_output real(DP), dimension(collision_system%fragments%nbody) :: vnoise type(lambda_obj) :: objective_function character(len=STRMAX) :: message @@ -584,7 +584,8 @@ subroutine fraggle_generate_rad_vel(collision_system, lfailure) objective_function = lambda_obj(radial_objective_function) tol = TOL_INIT do while(tol < TOL_MIN) - call util_minimize_bfgs(objective_function, nfrag, v_r_initial, tol, MAXLOOP, lfailure, fragments%v_r_mag) + call swiftest_util_minimize_bfgs(objective_function, nfrag, v_r_initial, tol, MAXLOOP, lfailure, v_r_output) + fragments%v_r_mag(1:nfrag) = v_r_output(1:nfrag) if (.not.lfailure) exit tol = tol * 2 ! Keep increasing the tolerance until we converge on a solution v_r_initial(:) = fragments%v_r_mag(:) diff --git a/src/main/swiftest_driver.f90 b/src/main/swiftest_driver.f90 index b7a0b330b..2b2229a3a 100644 --- a/src/main/swiftest_driver.f90 +++ b/src/main/swiftest_driver.f90 @@ -170,5 +170,5 @@ program swiftest_driver end associate end associate - call util_exit(SUCCESS) + call swiftest_util_exit(SUCCESS) end program swiftest_driver diff --git a/src/modules/base.f90 b/src/modules/base.f90 index 83e26eef6..ea28a5e18 100644 --- a/src/modules/base.f90 +++ b/src/modules/base.f90 @@ -387,7 +387,7 @@ subroutine netcdf_check(status, call_identifier) if(status /= nf90_noerr) then if (present(call_identifier)) write(*,*) "NetCDF error in ",trim(call_identifier) write(*,*) trim(nf90_strerror(status)) - call util_exit(FAILURE) + call swiftest_util_exit(FAILURE) end if return diff --git a/src/modules/collision.f90 b/src/modules/collision.f90 index bee281fd3..a4e6cd9bb 100644 --- a/src/modules/collision.f90 +++ b/src/modules/collision.f90 @@ -7,6 +7,7 @@ !! You should have received a copy of the GNU General Public License along with Swiftest. !! If not, see: https://www.gnu.org/licenses. + module collision !! author: The Purdue Swiftest Team - David A. Minton, Carlisle A. Wishard, Jennifer L.L. Pouplin, and Jacob R. Elliott !! diff --git a/src/modules/encounter.f90 b/src/modules/encounter.f90 index d47fca9ee..04e7170c6 100644 --- a/src/modules/encounter.f90 +++ b/src/modules/encounter.f90 @@ -22,7 +22,7 @@ module encounter integer(I8B) :: nenc = 0 !! Total number of encounters real(DP) :: t !! Time of encounter logical :: lcollision !! Indicates if the encounter resulted in at least one collision - real(DP), dimension(:), allocatable :: tcollision!! Time of collision + real(DP), dimension(:), allocatable :: tcollision !! Time of collision logical, dimension(:), allocatable :: lclosest !! indicates that thie pair of bodies is in currently at its closest approach point logical, dimension(:), allocatable :: lvdotr !! relative vdotr flag integer(I4B), dimension(:), allocatable :: status !! status of the interaction @@ -34,6 +34,7 @@ module encounter real(DP), dimension(:,:), allocatable :: r2 !! the position of body 2 in the encounter real(DP), dimension(:,:), allocatable :: v1 !! the velocity of body 1 in the encounter real(DP), dimension(:,:), allocatable :: v2 !! the velocity of body 2 in the encounter + integer(I4B), dimension(:), allocatable :: level !! Recursion level (used in SyMBA) contains procedure :: setup => encounter_setup_list !! A constructor that sets the number of encounters and allocates and initializes all arrays procedure :: append => encounter_util_append_list !! Appends elements from one structure to another diff --git a/src/modules/swiftest.f90 b/src/modules/swiftest.f90 index bd329f0d9..01afd9134 100644 --- a/src/modules/swiftest.f90 +++ b/src/modules/swiftest.f90 @@ -1060,7 +1060,7 @@ module subroutine swiftest_user_kick_getacch_body(self, system, param, t, lbeg) end subroutine swiftest_user_kick_getacch_body end interface - interface util_append + interface swiftest_util_append module subroutine swiftest_util_append_arr_char_string(arr, source, nold, nsrc, lsource_mask) implicit none character(len=STRMAX), dimension(:), allocatable, intent(inout) :: arr !! Destination array @@ -1272,7 +1272,7 @@ module subroutine swiftest_util_fill_tp(self, inserts, lfill_list) end subroutine swiftest_util_fill_tp end interface - interface util_fill + interface swiftest_util_fill module subroutine swiftest_util_fill_arr_char_string(keeps, inserts, lfill_list) implicit none character(len=STRMAX), dimension(:), allocatable, intent(inout) :: keeps !! Array of values to keep @@ -1430,7 +1430,7 @@ end subroutine swiftest_util_reset_storage end interface - interface util_resize + interface swiftest_util_resize module subroutine swiftest_util_resize_arr_char_string(arr, nnew) implicit none character(len=STRMAX), dimension(:), allocatable, intent(inout) :: arr !! Array to resize @@ -1588,28 +1588,28 @@ end subroutine swiftest_util_snapshot_system end interface - interface util_solve_linear_system - module function util_solve_linear_system_d(A,b,n,lerr) result(x) + interface swiftest_util_solve_linear_system + module function swiftest_util_solve_linear_system_d(A,b,n,lerr) result(x) implicit none integer(I4B), intent(in) :: n real(DP), dimension(:,:), intent(in) :: A real(DP), dimension(:), intent(in) :: b logical, intent(out) :: lerr real(DP), dimension(n) :: x - end function util_solve_linear_system_d + end function swiftest_util_solve_linear_system_d - module function util_solve_linear_system_q(A,b,n,lerr) result(x) + module function swiftest_util_solve_linear_system_q(A,b,n,lerr) result(x) implicit none integer(I4B), intent(in) :: n real(QP), dimension(:,:), intent(in) :: A real(QP), dimension(:), intent(in) :: b logical, intent(out) :: lerr real(QP), dimension(n) :: x - end function util_solve_linear_system_q + end function swiftest_util_solve_linear_system_q end interface interface - module function util_solve_rkf45(f, y0in, t1, dt0, tol) result(y1) + module function swiftest_util_solve_rkf45(f, y0in, t1, dt0, tol) result(y1) use lambda_function implicit none class(lambda_obj), intent(inout) :: f !! lambda function object that has been initialized to be a function of derivatives. The object will return with components lastarg and lasteval set @@ -1618,10 +1618,10 @@ module function util_solve_rkf45(f, y0in, t1, dt0, tol) result(y1) real(DP), intent(in) :: dt0 !! Initial step size guess real(DP), intent(in) :: tol !! Tolerance on solution real(DP), dimension(:), allocatable :: y1 !! Final result - end function util_solve_rkf45 + end function swiftest_util_solve_rkf45 end interface - interface util_sort + interface swiftest_util_sort pure module subroutine swiftest_util_sort_i4b(arr) implicit none integer(I4B), dimension(:), intent(inout) :: arr @@ -1666,9 +1666,9 @@ pure module subroutine swiftest_util_sort_index_dp(arr,ind) real(DP), dimension(:), intent(in) :: arr integer(I4B), dimension(:), allocatable, intent(inout) :: ind end subroutine swiftest_util_sort_index_dp - end interface util_sort + end interface swiftest_util_sort - interface util_sort_rearrange + interface swiftest_util_sort_rearrange pure module subroutine swiftest_util_sort_rearrange_arr_char_string(arr, ind, n) implicit none character(len=STRMAX), dimension(:), allocatable, intent(inout) :: arr !! Destination array @@ -1731,7 +1731,7 @@ pure module subroutine swiftest_util_sort_rearrange_arr_logical_I8Bind(arr, ind, integer(I8B), dimension(:), intent(in) :: ind !! Index to rearrange against integer(I8B), intent(in) :: n !! Number of elements in arr and ind to rearrange end subroutine swiftest_util_sort_rearrange_arr_logical_I8Bind - end interface util_sort_rearrange + end interface swiftest_util_sort_rearrange interface module subroutine swiftest_util_sort_rearrange_body(self, ind) @@ -1775,7 +1775,7 @@ end subroutine swiftest_util_sort_tp end interface - interface util_spill + interface swiftest_util_spill module subroutine swiftest_util_spill_arr_char_string(keeps, discards, lspill_list, ldestructive) implicit none character(len=STRMAX), dimension(:), allocatable, intent(inout) :: keeps !! Array of values to keep @@ -1868,7 +1868,7 @@ end subroutine swiftest_util_spill_tp end interface - interface util_unique + interface swiftest_util_unique module subroutine swiftest_util_unique_DP(input_array, output_array, index_map) implicit none real(DP), dimension(:), intent(in) :: input_array !! Unsorted input array @@ -1882,7 +1882,7 @@ module subroutine swiftest_util_unique_I4B(input_array, output_array, index_map) integer(I4B), dimension(:), allocatable, intent(out) :: output_array !! Sorted array of unique values integer(I4B), dimension(:), allocatable, intent(out) :: index_map !! An array of the same size as input_array that such that any for any index i, output_array(index_map(i)) = input_array(i) end subroutine swiftest_util_unique_I4B - end interface util_unique + end interface swiftest_util_unique interface module subroutine swiftest_util_valid_id_system(self, param) diff --git a/src/modules/symba.f90 b/src/modules/symba.f90 index 419a2b711..87359501c 100644 --- a/src/modules/symba.f90 +++ b/src/modules/symba.f90 @@ -76,18 +76,21 @@ module symba !> SyMBA class for tracking close encounters in a step - type, extends(encounter_list) :: symba_encounter - integer(I4B), dimension(:), allocatable :: level !! encounter recursion level + type, extends(collision_list_plpl) :: symba_list_plpl contains - procedure :: encounter_check => symba_encounter_check !! Checks if massive bodies are going through close encounters with each other - procedure :: kick => symba_kick_encounter !! Kick barycentric velocities of active test particles within SyMBA recursion - procedure :: setup => symba_setup_encounter_list !! A constructor that sets the number of encounters and allocates and initializes all arrays - procedure :: copy => symba_util_copy_encounter_list !! Copies elements from the source encounter list into self. - procedure :: dealloc => symba_util_dealloc_encounter_list !! Deallocates all allocatable arrays - procedure :: spill => symba_util_spill_encounter_list !! "Spills" bodies from one object to another depending on the results of a mask (uses the PACK intrinsic) - procedure :: append => symba_util_append_encounter_list !! Appends elements from one structure to another - final :: symba_util_final_encounter_list !! Finalizes the SyMBA test particle object - deallocates all allocatables - end type symba_encounter + procedure :: encounter_check => symba_encounter_check_list_plpl !! Checks if massive bodies are going through close encounters with each other + procedure :: kick => symba_kick_list_plpl !! Kick barycentric velocities of active massive bodies within SyMBA recursion + final :: symba_util_final_list_plpl !! Finalizes the SyMBA test particle object - deallocates all allocatables + end type symba_list_plpl + + + !> SyMBA class for tracking close encounters in a step + type, extends(collision_list_pltp) :: symba_list_pltp + contains + procedure :: encounter_check => symba_encounter_check_list_pltp !! Checks if massive bodies are going through close encounters with test particles + procedure :: kick => symba_kick_list_pltp !! Kick barycentric velocities of active test particles within SyMBA recursion + final :: symba_util_final_list_pltp !! Finalizes the SyMBA test particle object - deallocates all allocatables + end type symba_list_pltp type, extends(helio_nbody_system) :: symba_nbody_system @@ -102,7 +105,6 @@ module symba final :: symba_util_final_system !! Finalizes the SyMBA nbody system object - deallocates all allocatables end type symba_nbody_system - interface module subroutine symba_discard_pl(self, system, param) implicit none @@ -137,15 +139,25 @@ module function symba_encounter_check_pl(self, param, system, dt, irec) result(l logical :: lany_encounter !! Returns true if there is at least one close encounter end function symba_encounter_check_pl - module function symba_encounter_check(self, param, system, dt, irec) result(lany_encounter) + module function symba_encounter_check_list_plpl(self, param, system, dt, irec) result(lany_encounter) + implicit none + class(symba_list_plpl), intent(inout) :: self !! SyMBA pl-pl encounter list object + class(swiftest_parameters), intent(inout) :: param !! Current swiftest run configuration parameters + class(symba_nbody_system), intent(inout) :: system !! SyMBA nbody system object + real(DP), intent(in) :: dt !! step size + integer(I4B), intent(in) :: irec !! Current recursion level + logical :: lany_encounter !! Returns true if there is at least one close encounter + end function symba_encounter_check_list_plpl + + module function symba_encounter_check_list_pltp(self, param, system, dt, irec) result(lany_encounter) implicit none - class(symba_encounter), intent(inout) :: self !! SyMBA pl-pl encounter list object + class(symba_list_pltp), intent(inout) :: self !! SyMBA pl-tp encounter list object class(swiftest_parameters), intent(inout) :: param !! Current swiftest run configuration parameters class(symba_nbody_system), intent(inout) :: system !! SyMBA nbody system object real(DP), intent(in) :: dt !! step size integer(I4B), intent(in) :: irec !! Current recursion level logical :: lany_encounter !! Returns true if there is at least one close encounter - end function symba_encounter_check + end function symba_encounter_check_list_pltp module function symba_encounter_check_tp(self, param, system, dt, irec) result(lany_encounter) implicit none @@ -214,14 +226,23 @@ module subroutine symba_kick_getacch_tp(self, system, param, t, lbeg) logical, intent(in) :: lbeg !! Logical flag that determines whether or not this is the beginning or end of the step end subroutine symba_kick_getacch_tp - module subroutine symba_kick_encounter(self, system, dt, irec, sgn) + module subroutine symba_kick_list_plpl(self, system, dt, irec, sgn) implicit none - class(symba_encounter), intent(in) :: self !! SyMBA pl-tp encounter list object + class(symba_list_plpl), intent(in) :: self !! SyMBA pl-tp encounter list object class(symba_nbody_system), intent(inout) :: system !! SyMBA nbody system object real(DP), intent(in) :: dt !! step size integer(I4B), intent(in) :: irec !! Current recursion level integer(I4B), intent(in) :: sgn !! sign to be applied to acceleration - end subroutine symba_kick_encounter + end subroutine symba_kick_list_plpl + + module subroutine symba_kick_list_pltp(self, system, dt, irec, sgn) + implicit none + class(symba_list_pltp), intent(in) :: self !! SyMBA pl-tp encounter list object + class(symba_nbody_system), intent(inout) :: system !! SyMBA nbody system object + real(DP), intent(in) :: dt !! step size + integer(I4B), intent(in) :: irec !! Current recursion level + integer(I4B), intent(in) :: sgn !! sign to be applied to acceleration + end subroutine symba_kick_list_pltp module subroutine symba_setup_initialize_system(self, param) implicit none @@ -229,7 +250,6 @@ module subroutine symba_setup_initialize_system(self, param) class(swiftest_parameters), intent(inout) :: param !! Current run configuration parameters end subroutine symba_setup_initialize_system - module subroutine symba_setup_pl(self, n, param) implicit none class(symba_pl), intent(inout) :: self !! SyMBA massive body object @@ -237,12 +257,6 @@ module subroutine symba_setup_pl(self, n, param) class(swiftest_parameters), intent(in) :: param !! Current run configuration parameters end subroutine symba_setup_pl - module subroutine symba_setup_encounter_list(self,n) - implicit none - class(symba_encounter), intent(inout) :: self !! SyMBA pl-tp encounter structure - integer(I8B), intent(in) :: n !! Number of encounters to allocate space for - end subroutine symba_setup_encounter_list - module subroutine symba_setup_tp(self, n, param) implicit none class(symba_tp), intent(inout) :: self !! SyMBA test particle object @@ -286,13 +300,6 @@ module subroutine symba_step_reset_system(self, param) class(swiftest_parameters), intent(in) :: param !! Current run configuration parameters with SyMBA additions end subroutine symba_step_reset_system - module subroutine symba_util_append_encounter_list(self, source, lsource_mask) - implicit none - class(symba_encounter), intent(inout) :: self !! SyMBA encounter list object - class(encounter_list), intent(in) :: source !! Source object to append - logical, dimension(:), intent(in) :: lsource_mask !! Logical mask indicating which elements to append to - end subroutine symba_util_append_encounter_list - module subroutine symba_util_append_pl(self, source, lsource_mask) implicit none class(symba_pl), intent(inout) :: self !! SyMBA massive body object @@ -307,18 +314,6 @@ module subroutine symba_util_append_tp(self, source, lsource_mask) logical, dimension(:), intent(in) :: lsource_mask !! Logical mask indicating which elements to append to end subroutine symba_util_append_tp - module subroutine symba_util_copy_encounter_list(self, source) - use encounter, only : encounter_list - implicit none - class(symba_encounter), intent(inout) :: self !! Encounter list - class(encounter_list), intent(in) :: source !! Source object to copy into - end subroutine symba_util_copy_encounter_list - - module subroutine symba_util_dealloc_encounter_list(self) - implicit none - class(symba_encounter), intent(inout) :: self !! SyMBA encounter list - end subroutine symba_util_dealloc_encounter_list - module subroutine symba_util_dealloc_pl(self) implicit none class(symba_pl), intent(inout) :: self !! SyMBA massive body object @@ -352,10 +347,15 @@ module subroutine symba_util_flatten_eucl_plpl(self, param) class(swiftest_parameters), intent(inout) :: param !! Current run configuration parameters end subroutine symba_util_flatten_eucl_plpl - module subroutine symba_util_final_encounter_list(self) + module subroutine symba_util_final_list_plpl(self) implicit none - type(symba_encounter), intent(inout) :: self !! SyMBA encounter list object - end subroutine symba_util_final_encounter_list + type(symba_list_plpl), intent(inout) :: self !! SyMBA encounter list object + end subroutine symba_util_final_list_plpl + + module subroutine symba_util_final_list_pltp(self) + implicit none + type(symba_list_pltp), intent(inout) :: self !! SyMBA encounter list object + end subroutine symba_util_final_list_pltp module subroutine symba_util_final_pl(self) implicit none @@ -432,15 +432,6 @@ module subroutine symba_util_spill_pl(self, discards, lspill_list, ldestructive) logical, intent(in) :: ldestructive !! Logical flag indicating whether or not this operation should alter the keeps array or not end subroutine symba_util_spill_pl - module subroutine symba_util_spill_encounter_list(self, discards, lspill_list, ldestructive) - use encounter, only : encounter_list - implicit none - class(symba_encounter), intent(inout) :: self !! SyMBA pl-tp encounter list - class(encounter_list), intent(inout) :: discards !! Discarded object - logical, dimension(:), intent(in) :: lspill_list !! Logical array of bodies to spill into the discards - logical, intent(in) :: ldestructive !! Logical flag indicating whether or not this operation should alter body by removing the discard list - end subroutine symba_util_spill_encounter_list - module subroutine symba_util_spill_tp(self, discards, lspill_list, ldestructive) implicit none class(symba_tp), intent(inout) :: self !! SyMBA test particle object diff --git a/src/rmvs/rmvs_step.f90 b/src/rmvs/rmvs_step.f90 index 7c602a20b..c6f4bbab1 100644 --- a/src/rmvs/rmvs_step.f90 +++ b/src/rmvs/rmvs_step.f90 @@ -117,7 +117,7 @@ subroutine rmvs_interp_out(cb, pl, dt) write(*, *) xtmp(:,i) write(*, *) vtmp(:,i) write(*, *) " STOPPING " - call util_exit(FAILURE) + call swiftest_util_exit(FAILURE) end if end do end if @@ -139,7 +139,7 @@ subroutine rmvs_interp_out(cb, pl, dt) write(*, *) xtmp(:,i) write(*, *) vtmp(:,i) write(*, *) " STOPPING " - call util_exit(FAILURE) + call swiftest_util_exit(FAILURE) end if end do end if @@ -284,7 +284,7 @@ subroutine rmvs_interp_in(cb, pl, system, param, dt, outer_index) write(*, *) xtmp(:,i) write(*, *) vtmp(:,i) write(*, *) " STOPPING " - call util_exit(failure) + call swiftest_util_exit(failure) end if end do end if @@ -308,7 +308,7 @@ subroutine rmvs_interp_in(cb, pl, system, param, dt, outer_index) write(*, *) xtmp(:,i) write(*, *) vtmp(:,i) write(*, *) " STOPPING " - call util_exit(failure) + call swiftest_util_exit(failure) end if end do end if diff --git a/src/rmvs/rmvs_util.f90 b/src/rmvs/rmvs_util.f90 index a844507dc..09159afb5 100644 --- a/src/rmvs/rmvs_util.f90 +++ b/src/rmvs/rmvs_util.f90 @@ -25,21 +25,21 @@ module subroutine rmvs_util_append_pl(self, source, lsource_mask) select type(source) class is (rmvs_pl) associate(nold => self%nbody, nsrc => source%nbody) - call util_append(self%nenc, source%nenc, nold, nsrc, lsource_mask) - call util_append(self%tpenc1P, source%tpenc1P, nold, nsrc, lsource_mask) - call util_append(self%plind, source%plind, nold, nsrc, lsource_mask) + call swiftest_util_append(self%nenc, source%nenc, nold, nsrc, lsource_mask) + call swiftest_util_append(self%tpenc1P, source%tpenc1P, nold, nsrc, lsource_mask) + call swiftest_util_append(self%plind, source%plind, nold, nsrc, lsource_mask) ! The following are not implemented as RMVS doesn't make use of fill operations on pl type ! So they are here as a placeholder in case someone wants to extend the RMVS class for some reason - !call util_append(self%outer, source%outer, nold, nsrc, lsource_mask) - !call util_append(self%inner, source%inner, nold, nsrc, lsource_mask) - !call util_append(self%planetocentric, source%planetocentric, nold, nsrc, lsource_mask) + !call swiftest_util_append(self%outer, source%outer, nold, nsrc, lsource_mask) + !call swiftest_util_append(self%inner, source%inner, nold, nsrc, lsource_mask) + !call swiftest_util_append(self%planetocentric, source%planetocentric, nold, nsrc, lsource_mask) call whm_util_append_pl(self, source, lsource_mask) end associate class default write(*,*) "Invalid object passed to the append method. Source must be of class rmvs_pl or its descendents!" - call util_exit(FAILURE) + call swiftest_util_exit(FAILURE) end select return @@ -60,15 +60,15 @@ module subroutine rmvs_util_append_tp(self, source, lsource_mask) select type(source) class is (rmvs_tp) associate(nold => self%nbody, nsrc => source%nbody) - call util_append(self%lperi, source%lperi, nold, nsrc, lsource_mask) - call util_append(self%plperP, source%plperP, nold, nsrc, lsource_mask) - call util_append(self%plencP, source%plencP, nold, nsrc, lsource_mask) + call swiftest_util_append(self%lperi, source%lperi, nold, nsrc, lsource_mask) + call swiftest_util_append(self%plperP, source%plperP, nold, nsrc, lsource_mask) + call swiftest_util_append(self%plencP, source%plencP, nold, nsrc, lsource_mask) - call util_append_tp(self, source, lsource_mask) ! Note: whm_tp does not have its own append method, so we skip back to the base class + call swiftest_util_append_tp(self, source, lsource_mask) ! Note: whm_tp does not have its own append method, so we skip back to the base class end associate class default write(*,*) "Invalid object passed to the append method. Source must be of class rmvs_tp or its descendents!" - call util_exit(FAILURE) + call swiftest_util_exit(FAILURE) end select return @@ -140,7 +140,7 @@ module subroutine rmvs_util_dealloc_tp(self) if (allocated(self%rheliocentric)) deallocate(self%rheliocentric) call self%cb_heliocentric%dealloc() - call util_dealloc_tp(self) + call swiftest_util_dealloc_tp(self) return end subroutine rmvs_util_dealloc_tp @@ -161,20 +161,20 @@ module subroutine rmvs_util_fill_pl(self, inserts, lfill_list) associate(keeps => self) select type(inserts) class is (rmvs_pl) - call util_fill(keeps%nenc, inserts%nenc, lfill_list) - call util_fill(keeps%tpenc1P, inserts%tpenc1P, lfill_list) - call util_fill(keeps%plind, inserts%plind, lfill_list) + call swiftest_util_fill(keeps%nenc, inserts%nenc, lfill_list) + call swiftest_util_fill(keeps%tpenc1P, inserts%tpenc1P, lfill_list) + call swiftest_util_fill(keeps%plind, inserts%plind, lfill_list) ! The following are not implemented as RMVS doesn't make use of fill operations on pl type ! So they are here as a placeholder in case someone wants to extend the RMVS class for some reason - !call util_fill(keeps%outer, inserts%outer, lfill_list) - !call util_fill(keeps%inner, inserts%inner, lfill_list) - !call util_fill(keeps%planetocentric, inserts%planetocentric, lfill_list) + !call swiftest_util_fill(keeps%outer, inserts%outer, lfill_list) + !call swiftest_util_fill(keeps%inner, inserts%inner, lfill_list) + !call swiftest_util_fill(keeps%planetocentric, inserts%planetocentric, lfill_list) call whm_util_fill_pl(keeps, inserts, lfill_list) class default write(*,*) "Invalid object passed to the fill method. Source must be of class rmvs_pl or its descendents!" - call util_exit(FAILURE) + call swiftest_util_exit(FAILURE) end select end associate @@ -268,14 +268,14 @@ module subroutine rmvs_util_fill_tp(self, inserts, lfill_list) associate(keeps => self) select type(inserts) class is (rmvs_tp) - call util_fill(keeps%lperi, inserts%lperi, lfill_list) - call util_fill(keeps%plperP, inserts%plperP, lfill_list) - call util_fill(keeps%plencP, inserts%plencP, lfill_list) + call swiftest_util_fill(keeps%lperi, inserts%lperi, lfill_list) + call swiftest_util_fill(keeps%plperP, inserts%plperP, lfill_list) + call swiftest_util_fill(keeps%plencP, inserts%plencP, lfill_list) - call util_fill_tp(keeps, inserts, lfill_list) ! Note: whm_tp does not have its own fill method, so we skip back to the base class + call swiftest_util_fill_tp(keeps, inserts, lfill_list) ! Note: whm_tp does not have its own fill method, so we skip back to the base class class default write(*,*) "Invalid object passed to the fill method. Source must be of class rmvs_tp or its descendents!" - call util_exit(FAILURE) + call swiftest_util_exit(FAILURE) end select end associate @@ -292,15 +292,15 @@ module subroutine rmvs_util_resize_pl(self, nnew) class(rmvs_pl), intent(inout) :: self !! RMVS massive body object integer(I4B), intent(in) :: nnew !! New size neded - call util_resize(self%nenc, nnew) - call util_resize(self%tpenc1P, nnew) - call util_resize(self%plind, nnew) + call swiftest_util_resize(self%nenc, nnew) + call swiftest_util_resize(self%tpenc1P, nnew) + call swiftest_util_resize(self%plind, nnew) ! The following are not implemented as RMVS doesn't make use of resize operations on pl type ! So they are here as a placeholder in case someone wants to extend the RMVS class for some reason - !call util_resize(self%outer, nnew) - !call util_resize(self%inner, nnew) - !call util_resize(self%planetocentric, nnew) + !call swiftest_util_resize(self%outer, nnew) + !call swiftest_util_resize(self%inner, nnew) + !call swiftest_util_resize(self%planetocentric, nnew) call whm_util_resize_pl(self, nnew) return @@ -316,12 +316,12 @@ module subroutine rmvs_util_resize_tp(self, nnew) class(rmvs_tp), intent(inout) :: self !! RMVS test particle object integer(I4B), intent(in) :: nnew !! New size neded - call util_resize(self%lperi, nnew) - call util_resize(self%plperP, nnew) - call util_resize(self%plencP, nnew) - call util_resize(self%rheliocentric, nnew) + call swiftest_util_resize(self%lperi, nnew) + call swiftest_util_resize(self%plperP, nnew) + call swiftest_util_resize(self%plencP, nnew) + call swiftest_util_resize(self%rheliocentric, nnew) - call util_resize_tp(self, nnew) + call swiftest_util_resize_tp(self, nnew) return end subroutine rmvs_util_resize_tp @@ -352,11 +352,11 @@ module subroutine rmvs_util_sort_pl(self, sortby, ascending) associate(pl => self, npl => self%nbody) select case(sortby) case("nenc") - call util_sort(direction * pl%nenc(1:npl), ind) + call swiftest_util_sort(direction * pl%nenc(1:npl), ind) case("tpenc1P") - call util_sort(direction * pl%tpenc1P(1:npl), ind) + call swiftest_util_sort(direction * pl%tpenc1P(1:npl), ind) case("plind") - call util_sort(direction * pl%plind(1:npl), ind) + call swiftest_util_sort(direction * pl%plind(1:npl), ind) case("outer", "inner", "planetocentric", "lplanetocentric") write(*,*) 'Cannot sort by ' // trim(adjustl(sortby)) // '. Component not sortable!' case default ! Look for components in the parent class @@ -396,13 +396,13 @@ module subroutine rmvs_util_sort_tp(self, sortby, ascending) associate(tp => self, ntp => self%nbody) select case(sortby) case("plperP") - call util_sort(direction * tp%plperP(1:ntp), ind) + call swiftest_util_sort(direction * tp%plperP(1:ntp), ind) case("plencP") - call util_sort(direction * tp%plencP(1:ntp), ind) + call swiftest_util_sort(direction * tp%plencP(1:ntp), ind) case("lperi", "cb_heliocentric", "rheliocentric", "index", "ipleP", "lplanetocentric") write(*,*) 'Cannot sort by ' // trim(adjustl(sortby)) // '. Component not sortable!' case default ! Look for components in the parent class (*NOTE whm_tp does not need its own sort method, so we go straight to the swiftest_tp method) - call util_sort_tp(tp, sortby, ascending) + call swiftest_util_sort_tp(tp, sortby, ascending) return end select @@ -425,10 +425,10 @@ module subroutine rmvs_util_sort_rearrange_pl(self, ind) if (self%nbody == 0) return associate(pl => self, npl => self%nbody) - call util_sort_rearrange(pl%nenc, ind, npl) - call util_sort_rearrange(pl%tpenc1P, ind, npl) - call util_sort_rearrange(pl%plind, ind, npl) - call util_sort_rearrange_pl(pl,ind) + call swiftest_util_sort_rearrange(pl%nenc, ind, npl) + call swiftest_util_sort_rearrange(pl%tpenc1P, ind, npl) + call swiftest_util_sort_rearrange(pl%plind, ind, npl) + call swiftest_util_sort_rearrange_pl(pl,ind) end associate return @@ -448,11 +448,11 @@ module subroutine rmvs_util_sort_rearrange_tp(self, ind) if (self%nbody == 0) return associate(tp => self, ntp => self%nbody) - call util_sort_rearrange(tp%lperi, ind, ntp) - call util_sort_rearrange(tp%plperP, ind, ntp) - call util_sort_rearrange(tp%plencP, ind, ntp) - call util_sort_rearrange(tp%rheliocentric, ind, ntp) - call util_sort_rearrange_tp(tp,ind) + call swiftest_util_sort_rearrange(tp%lperi, ind, ntp) + call swiftest_util_sort_rearrange(tp%plperP, ind, ntp) + call swiftest_util_sort_rearrange(tp%plencP, ind, ntp) + call swiftest_util_sort_rearrange(tp%rheliocentric, ind, ntp) + call swiftest_util_sort_rearrange_tp(tp,ind) end associate return @@ -475,14 +475,14 @@ module subroutine rmvs_util_spill_pl(self, discards, lspill_list, ldestructive) associate(keeps => self) select type(discards) class is (rmvs_pl) - call util_spill(keeps%nenc, discards%nenc, lspill_list, ldestructive) - call util_spill(keeps%tpenc1P, discards%tpenc1P, lspill_list, ldestructive) - call util_spill(keeps%plind, discards%plind, lspill_list, ldestructive) + call swiftest_util_spill(keeps%nenc, discards%nenc, lspill_list, ldestructive) + call swiftest_util_spill(keeps%tpenc1P, discards%tpenc1P, lspill_list, ldestructive) + call swiftest_util_spill(keeps%plind, discards%plind, lspill_list, ldestructive) call whm_util_spill_pl(keeps, discards, lspill_list, ldestructive) class default write(*,*) "Invalid object passed to the spill method. Source must be of class rmvs_pl or its descendents!" - call util_exit(FAILURE) + call swiftest_util_exit(FAILURE) end select end associate @@ -506,14 +506,14 @@ module subroutine rmvs_util_spill_tp(self, discards, lspill_list, ldestructive) associate(keeps => self) select type(discards) class is (rmvs_tp) - call util_spill(keeps%lperi, discards%lperi, lspill_list, ldestructive) - call util_spill(keeps%plperP, discards%plperP, lspill_list, ldestructive) - call util_spill(keeps%plencP, discards%plencP, lspill_list, ldestructive) + call swiftest_util_spill(keeps%lperi, discards%lperi, lspill_list, ldestructive) + call swiftest_util_spill(keeps%plperP, discards%plperP, lspill_list, ldestructive) + call swiftest_util_spill(keeps%plencP, discards%plencP, lspill_list, ldestructive) - call util_spill_tp(keeps, discards, lspill_list, ldestructive) + call swiftest_util_spill_tp(keeps, discards, lspill_list, ldestructive) class default write(*,*) "Invalid object passed to the spill method. Source must be of class rmvs_tp or its descendents!" - call util_exit(FAILURE) + call swiftest_util_exit(FAILURE) end select end associate diff --git a/src/swiftest_procedures/swiftest_io.f90 b/src/swiftest_procedures/swiftest_io.f90 index 885e77985..581e7e001 100644 --- a/src/swiftest_procedures/swiftest_io.f90 +++ b/src/swiftest_procedures/swiftest_io.f90 @@ -178,7 +178,7 @@ module subroutine swiftest_io_conservation_report(self, param, lterminal) param%ioutput = param%ioutput + 1 call self%write_frame(nc, param) call nc%close() - call util_exit(FAILURE) + call swiftest_util_exit(FAILURE) end if end if end associate @@ -187,7 +187,7 @@ module subroutine swiftest_io_conservation_report(self, param, lterminal) 667 continue write(*,*) "Error writing energy and momentum tracking file: " // trim(adjustl(errmsg)) - call util_exit(FAILURE) + call swiftest_util_exit(FAILURE) end subroutine swiftest_io_conservation_report @@ -218,7 +218,7 @@ module subroutine swiftest_io_dump_param(self, param_file_name) 667 continue write(*,*) "Error opening parameter dump file " // trim(adjustl(errmsg)) - call util_exit(FAILURE) + call swiftest_util_exit(FAILURE) end subroutine swiftest_io_dump_param @@ -331,18 +331,18 @@ module subroutine swiftest_io_get_args(integrator, param_file_name, display_styl do i = 1,narg call get_command_argument(i, arg(i), status = ierr(i)) end do - if (any(ierr /= 0)) call util_exit(USAGE) + if (any(ierr /= 0)) call swiftest_util_exit(USAGE) else - call util_exit(USAGE) + call swiftest_util_exit(USAGE) end if if (narg == 1) then if (arg(1) == '-v' .or. arg(1) == '--version') then - call util_version() + call swiftest_util_version() else if (arg(1) == '-h' .or. arg(1) == '--help') then - call util_exit(HELP) + call swiftest_util_exit(HELP) else - call util_exit(USAGE) + call swiftest_util_exit(USAGE) end if else if (narg >= 2) then call io_toupper(arg(1)) @@ -366,7 +366,7 @@ module subroutine swiftest_io_get_args(integrator, param_file_name, display_styl case default integrator = UNKNOWN_INTEGRATOR write(*,*) trim(adjustl(arg(1))) // ' is not a valid integrator.' - call util_exit(USAGE) + call swiftest_util_exit(USAGE) end select param_file_name = trim(adjustl(arg(2))) end if @@ -377,7 +377,7 @@ module subroutine swiftest_io_get_args(integrator, param_file_name, display_styl call io_toupper(arg(3)) display_style = trim(adjustl(arg(3))) else - call util_exit(USAGE) + call swiftest_util_exit(USAGE) end if return @@ -1328,7 +1328,7 @@ module subroutine swiftest_io_read_in_cb(self, param) 667 continue write(*,*) "Error reading central body file: " // trim(adjustl(errmsg)) - call util_exit(FAILURE) + call swiftest_util_exit(FAILURE) end subroutine swiftest_io_read_in_cb @@ -1367,7 +1367,7 @@ module subroutine swiftest_io_read_in_system(self, param) end if ierr = self%read_frame(tmp_param%system_history%nc, tmp_param) deallocate(tmp_param) - if (ierr /=0) call util_exit(FAILURE) + if (ierr /=0) call swiftest_util_exit(FAILURE) end if param%loblatecb = ((self%cb%j2rp2 /= 0.0_DP) .or. (self%cb%j4rp4 /= 0.0_DP)) @@ -1481,7 +1481,7 @@ module function swiftest_io_read_frame_body(self, iu, param) result(ierr) class default write(*,*) "Error reading body file: " // trim(adjustl(errmsg)) end select - call util_exit(FAILURE) + call swiftest_util_exit(FAILURE) end function swiftest_io_read_frame_body @@ -1513,7 +1513,7 @@ module subroutine swiftest_io_read_in_param(self, param_file_name) 667 continue write(self%display_unit,*) "Error reading parameter file: " // trim(adjustl(errmsg)) - call util_exit(FAILURE) + call swiftest_util_exit(FAILURE) end subroutine swiftest_io_read_in_param @@ -1539,7 +1539,7 @@ module subroutine swiftest_io_set_display_param(self, display_style) self%log_output = .true. case default write(*,*) display_style, " is an unknown display style" - call util_exit(USAGE) + call swiftest_util_exit(USAGE) end select self%display_style = display_style @@ -1548,7 +1548,7 @@ module subroutine swiftest_io_set_display_param(self, display_style) 667 continue write(*,*) "Error opening swiftest log file: " // trim(adjustl(errmsg)) - call util_exit(FAILURE) + call swiftest_util_exit(FAILURE) end subroutine swiftest_io_set_display_param @@ -1650,7 +1650,7 @@ module subroutine swiftest_io_write_frame_system(self, param) 667 continue write(*,*) "Error writing system frame: " // trim(adjustl(errmsg)) - call util_exit(FAILURE) + call swiftest_util_exit(FAILURE) end subroutine swiftest_io_write_frame_system end submodule s_io diff --git a/src/swiftest_procedures/swiftest_io_netcdf.f90 b/src/swiftest_procedures/swiftest_io_netcdf.f90 index f45c29639..04f4b6b56 100644 --- a/src/swiftest_procedures/swiftest_io_netcdf.f90 +++ b/src/swiftest_procedures/swiftest_io_netcdf.f90 @@ -262,7 +262,7 @@ module subroutine swiftest_io_netcdf_initialize_output(self, param) 667 continue write(*,*) "Error creating NetCDF output file. " // trim(adjustl(errmsg)) - call util_exit(FAILURE) + call swiftest_util_exit(FAILURE) end subroutine swiftest_io_netcdf_initialize_output @@ -453,19 +453,19 @@ module function swiftest_io_netcdf_read_frame_system(self, nc, param) result(ier if (npl_check /= npl) then write(*,*) "Error reading in NetCDF file: The recorded value of npl does not match the number of active massive bodies" - call util_exit(failure) + call swiftest_util_exit(failure) end if if (ntp_check /= ntp) then write(*,*) "Error reading in NetCDF file: The recorded value of ntp does not match the number of active test particles" - call util_exit(failure) + call swiftest_util_exit(failure) end if if (param%integrator == INT_SYMBA) then nplm_check = count(pack(rtemp,plmask) > param%GMTINY ) if (nplm_check /= pl%nplm) then write(*,*) "Error reading in NetCDF file: The recorded value of nplm does not match the number of active fully interacting massive bodies" - call util_exit(failure) + call swiftest_util_exit(failure) end if end if @@ -980,7 +980,7 @@ module subroutine swiftest_io_netcdf_write_frame_body(self, nc, param) associate(n => self%nbody) if (n == 0) return - call util_sort(self%id(1:n), ind) + call swiftest_util_sort(self%id(1:n), ind) do i = 1, n j = ind(i) @@ -1166,7 +1166,7 @@ module subroutine swiftest_io_netcdf_write_info_body(self, nc, param) class is (swiftest_body) associate(n => self%nbody) if (n == 0) return - call util_sort(self%id(1:n), ind) + call swiftest_util_sort(self%id(1:n), ind) do i = 1, n j = ind(i) diff --git a/src/swiftest_procedures/swiftest_setup.f90 b/src/swiftest_procedures/swiftest_setup.f90 index 5e20ede68..0f3a274fb 100644 --- a/src/swiftest_procedures/swiftest_setup.f90 +++ b/src/swiftest_procedures/swiftest_setup.f90 @@ -111,7 +111,7 @@ module subroutine swiftest_setup_construct_system(system, param) write(*,*) 'RINGMOONS-SyMBA integrator not yet enabled' case default write(*,*) 'Unkown integrator',param%integrator - call util_exit(FAILURE) + call swiftest_util_exit(FAILURE) end select end select diff --git a/src/swiftest_procedures/swiftest_util.f90 b/src/swiftest_procedures/swiftest_util.f90 index 40f614f42..edc4c6693 100644 --- a/src/swiftest_procedures/swiftest_util.f90 +++ b/src/swiftest_procedures/swiftest_util.f90 @@ -29,7 +29,7 @@ module subroutine swiftest_util_append_arr_char_string(arr, source, nold, nsrc, if (.not.allocated(arr)) then allocate(arr(nold+nnew)) else - call util_resize(arr, nold + nnew) + call swiftest_util_resize(arr, nold + nnew) end if arr(nold + 1:nold + nnew) = pack(source(1:nsrc), lsource_mask(1:nsrc)) @@ -57,7 +57,7 @@ module subroutine swiftest_util_append_arr_DP(arr, source, nold, nsrc, lsource_m if (.not.allocated(arr)) then allocate(arr(nold+nnew)) else - call util_resize(arr, nold + nnew) + call swiftest_util_resize(arr, nold + nnew) end if arr(nold + 1:nold + nnew) = pack(source(1:nsrc), lsource_mask(1:nsrc)) @@ -85,7 +85,7 @@ module subroutine swiftest_util_append_arr_DPvec(arr, source, nold, nsrc, lsourc if (.not.allocated(arr)) then allocate(arr(NDIM,nold+nnew)) else - call util_resize(arr, nold + nnew) + call swiftest_util_resize(arr, nold + nnew) end if arr(1, nold + 1:nold + nnew) = pack(source(1,1:nsrc), lsource_mask(1:nsrc)) @@ -115,7 +115,7 @@ module subroutine swiftest_util_append_arr_I4B(arr, source, nold, nsrc, lsource_ if (.not.allocated(arr)) then allocate(arr(nold+nnew)) else - call util_resize(arr, nold + nnew) + call swiftest_util_resize(arr, nold + nnew) end if arr(nold + 1:nold + nnew) = pack(source(1:nsrc), lsource_mask(1:nsrc)) @@ -144,14 +144,14 @@ module subroutine swiftest_util_append_arr_info(arr, source, nold, nsrc, lsource if (.not.allocated(arr)) then allocate(arr(nold+nnew)) else - call util_resize(arr, nold + nnew) + call swiftest_util_resize(arr, nold + nnew) end if allocate(idx(nnew)) idx = pack([(i, i = 1, nsrc)], lsource_mask(1:nsrc)) - call util_copy_particle_info_arr(source(1:nsrc), arr(nold+1:nold+nnew), idx) + call swiftest_util_copy_particle_info_arr(source(1:nsrc), arr(nold+1:nold+nnew), idx) return end subroutine swiftest_util_append_arr_info @@ -176,7 +176,7 @@ module subroutine swiftest_util_append_arr_logical(arr, source, nold, nsrc, lsou if (.not.allocated(arr)) then allocate(arr(nold+nnew)) else - call util_resize(arr, nold + nnew) + call swiftest_util_resize(arr, nold + nnew) end if arr(nold + 1:nold + nnew) = pack(source(1:nsrc), lsource_mask(1:nsrc)) @@ -202,27 +202,27 @@ module subroutine swiftest_util_append_body(self, source, lsource_mask) nsrc = source%nbody nnew = count(lsource_mask(1:nsrc)) - call util_append(self%info, source%info, nold, nsrc, lsource_mask) - call util_append(self%id, source%id, nold, nsrc, lsource_mask) - call util_append(self%status, source%status, nold, nsrc, lsource_mask) - call util_append(self%ldiscard, source%ldiscard, nold, nsrc, lsource_mask) - call util_append(self%lmask, source%lmask, nold, nsrc, lsource_mask) - call util_append(self%mu, source%mu, nold, nsrc, lsource_mask) - call util_append(self%rh, source%rh, nold, nsrc, lsource_mask) - call util_append(self%vh, source%vh, nold, nsrc, lsource_mask) - call util_append(self%rb, source%rb, nold, nsrc, lsource_mask) - call util_append(self%vb, source%vb, nold, nsrc, lsource_mask) - call util_append(self%ah, source%ah, nold, nsrc, lsource_mask) - call util_append(self%aobl, source%aobl, nold, nsrc, lsource_mask) - call util_append(self%atide, source%atide, nold, nsrc, lsource_mask) - call util_append(self%agr, source%agr, nold, nsrc, lsource_mask) - call util_append(self%ir3h, source%ir3h, nold, nsrc, lsource_mask) - call util_append(self%a, source%a, nold, nsrc, lsource_mask) - call util_append(self%e, source%e, nold, nsrc, lsource_mask) - call util_append(self%inc, source%inc, nold, nsrc, lsource_mask) - call util_append(self%capom, source%capom, nold, nsrc, lsource_mask) - call util_append(self%omega, source%omega, nold, nsrc, lsource_mask) - call util_append(self%capm, source%capm, nold, nsrc, lsource_mask) + call swiftest_util_append(self%info, source%info, nold, nsrc, lsource_mask) + call swiftest_util_append(self%id, source%id, nold, nsrc, lsource_mask) + call swiftest_util_append(self%status, source%status, nold, nsrc, lsource_mask) + call swiftest_util_append(self%ldiscard, source%ldiscard, nold, nsrc, lsource_mask) + call swiftest_util_append(self%lmask, source%lmask, nold, nsrc, lsource_mask) + call swiftest_util_append(self%mu, source%mu, nold, nsrc, lsource_mask) + call swiftest_util_append(self%rh, source%rh, nold, nsrc, lsource_mask) + call swiftest_util_append(self%vh, source%vh, nold, nsrc, lsource_mask) + call swiftest_util_append(self%rb, source%rb, nold, nsrc, lsource_mask) + call swiftest_util_append(self%vb, source%vb, nold, nsrc, lsource_mask) + call swiftest_util_append(self%ah, source%ah, nold, nsrc, lsource_mask) + call swiftest_util_append(self%aobl, source%aobl, nold, nsrc, lsource_mask) + call swiftest_util_append(self%atide, source%atide, nold, nsrc, lsource_mask) + call swiftest_util_append(self%agr, source%agr, nold, nsrc, lsource_mask) + call swiftest_util_append(self%ir3h, source%ir3h, nold, nsrc, lsource_mask) + call swiftest_util_append(self%a, source%a, nold, nsrc, lsource_mask) + call swiftest_util_append(self%e, source%e, nold, nsrc, lsource_mask) + call swiftest_util_append(self%inc, source%inc, nold, nsrc, lsource_mask) + call swiftest_util_append(self%capom, source%capom, nold, nsrc, lsource_mask) + call swiftest_util_append(self%omega, source%omega, nold, nsrc, lsource_mask) + call swiftest_util_append(self%capm, source%capm, nold, nsrc, lsource_mask) self%nbody = nold + nnew @@ -244,28 +244,28 @@ module subroutine swiftest_util_append_pl(self, source, lsource_mask) select type(source) class is (swiftest_pl) associate(nold => self%nbody, nsrc => source%nbody) - call util_append(self%mass, source%mass, nold, nsrc, lsource_mask) - call util_append(self%Gmass, source%Gmass, nold, nsrc, lsource_mask) - call util_append(self%rhill, source%rhill, nold, nsrc, lsource_mask) - call util_append(self%renc, source%renc, nold, nsrc, lsource_mask) - call util_append(self%radius, source%radius, nold, nsrc, lsource_mask) - call util_append(self%rbeg, source%rbeg, nold, nsrc, lsource_mask) - call util_append(self%rend, source%rend, nold, nsrc, lsource_mask) - call util_append(self%vbeg, source%vbeg, nold, nsrc, lsource_mask) - call util_append(self%density, source%density, nold, nsrc, lsource_mask) - call util_append(self%Ip, source%Ip, nold, nsrc, lsource_mask) - call util_append(self%rot, source%rot, nold, nsrc, lsource_mask) - call util_append(self%k2, source%k2, nold, nsrc, lsource_mask) - call util_append(self%Q, source%Q, nold, nsrc, lsource_mask) - call util_append(self%tlag, source%tlag, nold, nsrc, lsource_mask) + call swiftest_util_append(self%mass, source%mass, nold, nsrc, lsource_mask) + call swiftest_util_append(self%Gmass, source%Gmass, nold, nsrc, lsource_mask) + call swiftest_util_append(self%rhill, source%rhill, nold, nsrc, lsource_mask) + call swiftest_util_append(self%renc, source%renc, nold, nsrc, lsource_mask) + call swiftest_util_append(self%radius, source%radius, nold, nsrc, lsource_mask) + call swiftest_util_append(self%rbeg, source%rbeg, nold, nsrc, lsource_mask) + call swiftest_util_append(self%rend, source%rend, nold, nsrc, lsource_mask) + call swiftest_util_append(self%vbeg, source%vbeg, nold, nsrc, lsource_mask) + call swiftest_util_append(self%density, source%density, nold, nsrc, lsource_mask) + call swiftest_util_append(self%Ip, source%Ip, nold, nsrc, lsource_mask) + call swiftest_util_append(self%rot, source%rot, nold, nsrc, lsource_mask) + call swiftest_util_append(self%k2, source%k2, nold, nsrc, lsource_mask) + call swiftest_util_append(self%Q, source%Q, nold, nsrc, lsource_mask) + call swiftest_util_append(self%tlag, source%tlag, nold, nsrc, lsource_mask) if (allocated(self%k_plpl)) deallocate(self%k_plpl) - call util_append_body(self, source, lsource_mask) + call swiftest_util_append_body(self, source, lsource_mask) end associate class default write(*,*) "Invalid object passed to the append method. Source must be of class swiftest_pl or its descendents" - call util_exit(FAILURE) + call swiftest_util_exit(FAILURE) end select return @@ -286,15 +286,15 @@ module subroutine swiftest_util_append_tp(self, source, lsource_mask) select type(source) class is (swiftest_tp) associate(nold => self%nbody, nsrc => source%nbody) - call util_append(self%isperi, source%isperi, nold, nsrc, lsource_mask) - call util_append(self%peri, source%peri, nold, nsrc, lsource_mask) - call util_append(self%atp, source%atp, nold, nsrc, lsource_mask) + call swiftest_util_append(self%isperi, source%isperi, nold, nsrc, lsource_mask) + call swiftest_util_append(self%peri, source%peri, nold, nsrc, lsource_mask) + call swiftest_util_append(self%atp, source%atp, nold, nsrc, lsource_mask) - call util_append_body(self, source, lsource_mask) + call swiftest_util_append_body(self, source, lsource_mask) end associate class default write(*,*) "Invalid object passed to the append method. Source must be of class swiftest_tp or its descendents" - call util_exit(FAILURE) + call swiftest_util_exit(FAILURE) end select return @@ -750,7 +750,7 @@ module subroutine swiftest_util_dealloc_pl(self) deallocate(self%kin) end if - call util_dealloc_body(self) + call swiftest_util_dealloc_body(self) return end subroutine swiftest_util_dealloc_pl @@ -770,7 +770,7 @@ module subroutine swiftest_util_dealloc_tp(self) if (allocated(self%atp)) deallocate(self%atp) if (allocated(self%k_pltp)) deallocate(self%k_pltp) - call util_dealloc_body(self) + call swiftest_util_dealloc_body(self) return end subroutine swiftest_util_dealloc_tp @@ -916,7 +916,7 @@ module subroutine swiftest_util_fill_arr_info(keeps, inserts, lfill_list) allocate(insert_idx(ninsert)) insert_idx(:) = pack([(i, i = 1, nkeep)], lfill_list) - call util_copy_particle_info_arr(inserts, keeps, insert_idx) + call swiftest_util_copy_particle_info_arr(inserts, keeps, insert_idx) return end subroutine swiftest_util_fill_arr_info @@ -956,32 +956,32 @@ module subroutine swiftest_util_fill_body(self, inserts, lfill_list) ! For each component, pack the discarded bodies into the discard object and do the inverse with the keeps !> Fill all the common components associate(keeps => self) - call util_fill(keeps%id, inserts%id, lfill_list) - call util_fill(keeps%info, inserts%info, lfill_list) - call util_fill(keeps%lmask, inserts%lmask, lfill_list) - call util_fill(keeps%status, inserts%status, lfill_list) - call util_fill(keeps%ldiscard, inserts%ldiscard, lfill_list) - call util_fill(keeps%lcollision, inserts%lcollision, lfill_list) - call util_fill(keeps%lencounter, inserts%lencounter, lfill_list) - call util_fill(keeps%mu, inserts%mu, lfill_list) - call util_fill(keeps%rh, inserts%rh, lfill_list) - call util_fill(keeps%vh, inserts%vh, lfill_list) - call util_fill(keeps%rb, inserts%rb, lfill_list) - call util_fill(keeps%vb, inserts%vb, lfill_list) - call util_fill(keeps%ah, inserts%ah, lfill_list) - call util_fill(keeps%aobl, inserts%aobl, lfill_list) - call util_fill(keeps%agr, inserts%agr, lfill_list) - call util_fill(keeps%atide, inserts%atide, lfill_list) - call util_fill(keeps%ir3h, inserts%ir3h, lfill_list) - call util_fill(keeps%isperi, inserts%isperi, lfill_list) - call util_fill(keeps%peri, inserts%peri, lfill_list) - call util_fill(keeps%atp, inserts%atp, lfill_list) - call util_fill(keeps%a, inserts%a, lfill_list) - call util_fill(keeps%e, inserts%e, lfill_list) - call util_fill(keeps%inc, inserts%inc, lfill_list) - call util_fill(keeps%capom, inserts%capom, lfill_list) - call util_fill(keeps%omega, inserts%omega, lfill_list) - call util_fill(keeps%capm, inserts%capm, lfill_list) + call swiftest_util_fill(keeps%id, inserts%id, lfill_list) + call swiftest_util_fill(keeps%info, inserts%info, lfill_list) + call swiftest_util_fill(keeps%lmask, inserts%lmask, lfill_list) + call swiftest_util_fill(keeps%status, inserts%status, lfill_list) + call swiftest_util_fill(keeps%ldiscard, inserts%ldiscard, lfill_list) + call swiftest_util_fill(keeps%lcollision, inserts%lcollision, lfill_list) + call swiftest_util_fill(keeps%lencounter, inserts%lencounter, lfill_list) + call swiftest_util_fill(keeps%mu, inserts%mu, lfill_list) + call swiftest_util_fill(keeps%rh, inserts%rh, lfill_list) + call swiftest_util_fill(keeps%vh, inserts%vh, lfill_list) + call swiftest_util_fill(keeps%rb, inserts%rb, lfill_list) + call swiftest_util_fill(keeps%vb, inserts%vb, lfill_list) + call swiftest_util_fill(keeps%ah, inserts%ah, lfill_list) + call swiftest_util_fill(keeps%aobl, inserts%aobl, lfill_list) + call swiftest_util_fill(keeps%agr, inserts%agr, lfill_list) + call swiftest_util_fill(keeps%atide, inserts%atide, lfill_list) + call swiftest_util_fill(keeps%ir3h, inserts%ir3h, lfill_list) + call swiftest_util_fill(keeps%isperi, inserts%isperi, lfill_list) + call swiftest_util_fill(keeps%peri, inserts%peri, lfill_list) + call swiftest_util_fill(keeps%atp, inserts%atp, lfill_list) + call swiftest_util_fill(keeps%a, inserts%a, lfill_list) + call swiftest_util_fill(keeps%e, inserts%e, lfill_list) + call swiftest_util_fill(keeps%inc, inserts%inc, lfill_list) + call swiftest_util_fill(keeps%capom, inserts%capom, lfill_list) + call swiftest_util_fill(keeps%omega, inserts%omega, lfill_list) + call swiftest_util_fill(keeps%capm, inserts%capm, lfill_list) ! This is the base class, so will be the last to be called in the cascade. keeps%nbody = size(keeps%id(:)) @@ -1007,27 +1007,27 @@ module subroutine swiftest_util_fill_pl(self, inserts, lfill_list) select type (inserts) ! The standard requires us to select the type of both arguments in order to access all the components class is (swiftest_pl) !> Fill components specific to the massive body class - call util_fill(keeps%mass, inserts%mass, lfill_list) - call util_fill(keeps%Gmass, inserts%Gmass, lfill_list) - call util_fill(keeps%rhill, inserts%rhill, lfill_list) - call util_fill(keeps%renc, inserts%renc, lfill_list) - call util_fill(keeps%radius, inserts%radius, lfill_list) - call util_fill(keeps%density, inserts%density, lfill_list) - call util_fill(keeps%rbeg, inserts%rbeg, lfill_list) - call util_fill(keeps%rend, inserts%rend, lfill_list) - call util_fill(keeps%vbeg, inserts%vbeg, lfill_list) - call util_fill(keeps%Ip, inserts%Ip, lfill_list) - call util_fill(keeps%rot, inserts%rot, lfill_list) - call util_fill(keeps%k2, inserts%k2, lfill_list) - call util_fill(keeps%Q, inserts%Q, lfill_list) - call util_fill(keeps%tlag, inserts%tlag, lfill_list) - call util_fill(keeps%kin, inserts%kin, lfill_list) - call util_fill(keeps%nplenc, inserts%nplenc, lfill_list) - call util_fill(keeps%ntpenc, inserts%ntpenc, lfill_list) + call swiftest_util_fill(keeps%mass, inserts%mass, lfill_list) + call swiftest_util_fill(keeps%Gmass, inserts%Gmass, lfill_list) + call swiftest_util_fill(keeps%rhill, inserts%rhill, lfill_list) + call swiftest_util_fill(keeps%renc, inserts%renc, lfill_list) + call swiftest_util_fill(keeps%radius, inserts%radius, lfill_list) + call swiftest_util_fill(keeps%density, inserts%density, lfill_list) + call swiftest_util_fill(keeps%rbeg, inserts%rbeg, lfill_list) + call swiftest_util_fill(keeps%rend, inserts%rend, lfill_list) + call swiftest_util_fill(keeps%vbeg, inserts%vbeg, lfill_list) + call swiftest_util_fill(keeps%Ip, inserts%Ip, lfill_list) + call swiftest_util_fill(keeps%rot, inserts%rot, lfill_list) + call swiftest_util_fill(keeps%k2, inserts%k2, lfill_list) + call swiftest_util_fill(keeps%Q, inserts%Q, lfill_list) + call swiftest_util_fill(keeps%tlag, inserts%tlag, lfill_list) + call swiftest_util_fill(keeps%kin, inserts%kin, lfill_list) + call swiftest_util_fill(keeps%nplenc, inserts%nplenc, lfill_list) + call swiftest_util_fill(keeps%ntpenc, inserts%ntpenc, lfill_list) if (allocated(keeps%k_plpl)) deallocate(keeps%k_plpl) - call util_fill_body(keeps, inserts, lfill_list) + call swiftest_util_fill_body(keeps, inserts, lfill_list) class default write(*,*) 'Error! fill method called for incompatible return type on swiftest_pl' end select @@ -1052,9 +1052,9 @@ module subroutine swiftest_util_fill_tp(self, inserts, lfill_list) select type(inserts) class is (swiftest_tp) !> Spill components specific to the test particle class - call util_fill(keeps%nplenc, inserts%nplenc, lfill_list) + call swiftest_util_fill(keeps%nplenc, inserts%nplenc, lfill_list) - call util_fill_body(keeps, inserts, lfill_list) + call swiftest_util_fill_body(keeps, inserts, lfill_list) class default write(*,*) 'Error! fill method called for incompatible return type on swiftest_tp' end select @@ -1186,7 +1186,7 @@ module subroutine swiftest_util_flatten_eucl_plpl(self, param) param%lflatten_interactions = .false. else do concurrent (i=1:npl, j=1:npl, j>i) - call util_flatten_eucl_ij_to_k(self%nbody, i, j, k) + call swiftest_util_flatten_eucl_ij_to_k(self%nbody, i, j, k) self%k_plpl(1, k) = i self%k_plpl(2, k) = j end do @@ -1312,9 +1312,9 @@ module subroutine swiftest_util_get_energy_momentum_system(self, param) end if if (param%lflatten_interactions) then - call util_get_energy_potential_flat(npl, pl%nplpl, pl%k_plpl, pl%lmask, cb%Gmass, pl%Gmass, pl%mass, pl%rb, system%pe) + call swiftest_util_get_energy_potential_flat(npl, pl%nplpl, pl%k_plpl, pl%lmask, cb%Gmass, pl%Gmass, pl%mass, pl%rb, system%pe) else - call util_get_energy_potential_triangular(npl, pl%lmask, cb%Gmass, pl%Gmass, pl%mass, pl%rb, system%pe) + call swiftest_util_get_energy_potential_triangular(npl, pl%lmask, cb%Gmass, pl%Gmass, pl%mass, pl%rb, system%pe) end if ! Potential energy from the oblateness term @@ -1577,12 +1577,12 @@ module subroutine swiftest_util_index_map_storage(self) integer(I4B), dimension(:), allocatable :: idvals real(DP), dimension(:), allocatable :: tvals - call util_get_vals_storage(self, idvals, tvals) + call swiftest_util_get_vals_storage(self, idvals, tvals) - call util_unique(idvals,self%idvals,self%idmap) + call swiftest_util_unique(idvals,self%idvals,self%idmap) self%nid = size(self%idvals) - call util_unique(tvals,self%tvals,self%tmap) + call swiftest_util_unique(tvals,self%tvals,self%tmap) self%nt = size(self%tvals) return @@ -2111,7 +2111,7 @@ subroutine swiftest_quadfit(f, x0, S, N, eps, lo, hi, lerr) lhs(2, :) = row_2 lhs(3, :) = row_3 ! Solve system of equations - soln(:) = util_solve_linear_system(lhs, rhs, 3, lerr) + soln(:) = swiftest_util_solve_linear_system(lhs, rhs, 3, lerr) call ieee_set_flag(ieee_all, .false.) ! Set all flags back to quiet call ieee_set_halting_mode(ieee_divide_by_zero, .false.) if (lerr) then @@ -2259,7 +2259,7 @@ module subroutine swiftest_util_rearray_pl(self, system, param) integer(I4B), dimension(:), allocatable :: levelg_orig_pl, levelm_orig_pl, levelg_orig_tp, levelm_orig_tp integer(I4B), dimension(:), allocatable :: nplenc_orig_pl, nplenc_orig_tp, ntpenc_orig_pl - associate(pl => self, pl_adds => system%pl_adds) + associate(pl => self, tp => system%tp, pl_adds => system%pl_adds) npl = pl%nbody nadd = pl_adds%nbody @@ -2327,32 +2327,38 @@ module subroutine swiftest_util_rearray_pl(self, system, param) ! Reset the kinship trackers call pl%reset_kinship([(i, i=1, npl)]) - ! Re-build the zero-level encounter list, being sure to save the original level information for all bodies - allocate(nplenc_orig_pl, source=pl%nplenc) + ! Re-build the encounter list + + + ! Be sure to get the level info if this is a SyMBA system + select type(system) + class is (symba_nbody_system) select type(pl) class is (symba_pl) + select type(tp) + class is (symba_tp) allocate(levelg_orig_pl, source=pl%levelg) allocate(levelm_orig_pl, source=pl%levelm) + allocate(nplenc_orig_tp, source=tp%nplenc) call move_alloc(levelg_orig_pl, pl%levelg) call move_alloc(levelm_orig_pl, pl%levelm) - end select - lencounter = pl%encounter_check(param, system, param%dt, system%irec) - if (system%tp%nbody > 0) then - allocate(ntpenc_orig_pl, source=pl%ntpenc) - allocate(nplenc_orig_tp, source=tp%nplenc) - select type(tp => system%tp) - class is (symba_tp) + call move_alloc(nplenc_orig_pl, pl%nplenc) + lencounter = pl%encounter_check(param, system, param%dt, system%irec) + if (tp%nbody > 0) then allocate(levelg_orig_tp, source=tp%levelg) allocate(levelm_orig_tp, source=tp%levelm) + allocate(nplenc_orig_tp, source=tp%nplenc) + allocate(ntpenc_orig_pl, source=pl%ntpenc) lencounter = tp%encounter_check(param, system, param%dt, system%irec) call move_alloc(levelg_orig_tp, tp%levelg) call move_alloc(levelm_orig_tp, tp%levelm) call move_alloc(nplenc_orig_tp, tp%nplenc) call move_alloc(ntpenc_orig_pl, pl%ntpenc) - end select - end if + end if + end select + end select + end select - call move_alloc(nplenc_orig_pl, pl%nplenc) ! Re-index the encounter list as the index values may have changed if (nenc_old > 0) then @@ -2689,9 +2695,9 @@ module subroutine swiftest_util_resize_arr_info(arr, nnew) allocate(tmp(nnew)) if (nnew > nold) then - call util_copy_particle_info_arr(arr(1:nold), tmp(1:nold)) + call swiftest_util_copy_particle_info_arr(arr(1:nold), tmp(1:nold)) else - call util_copy_particle_info_arr(arr(1:nnew), tmp(1:nnew)) + call swiftest_util_copy_particle_info_arr(arr(1:nnew), tmp(1:nnew)) end if call move_alloc(tmp, arr) @@ -2754,29 +2760,29 @@ module subroutine swiftest_util_resize_body(self, nnew) class(swiftest_body), intent(inout) :: self !! Swiftest body object integer(I4B), intent(in) :: nnew !! New size neded - call util_resize(self%info, nnew) - call util_resize(self%id, nnew) - call util_resize(self%status, nnew) - call util_resize(self%lcollision, nnew) - call util_resize(self%lencounter, nnew) - call util_resize(self%ldiscard, nnew) - call util_resize(self%lmask, nnew) - call util_resize(self%mu, nnew) - call util_resize(self%rh, nnew) - call util_resize(self%vh, nnew) - call util_resize(self%rb, nnew) - call util_resize(self%vb, nnew) - call util_resize(self%ah, nnew) - call util_resize(self%aobl, nnew) - call util_resize(self%atide, nnew) - call util_resize(self%agr, nnew) - call util_resize(self%ir3h, nnew) - call util_resize(self%a, nnew) - call util_resize(self%e, nnew) - call util_resize(self%inc, nnew) - call util_resize(self%capom, nnew) - call util_resize(self%omega, nnew) - call util_resize(self%capm, nnew) + call swiftest_util_resize(self%info, nnew) + call swiftest_util_resize(self%id, nnew) + call swiftest_util_resize(self%status, nnew) + call swiftest_util_resize(self%lcollision, nnew) + call swiftest_util_resize(self%lencounter, nnew) + call swiftest_util_resize(self%ldiscard, nnew) + call swiftest_util_resize(self%lmask, nnew) + call swiftest_util_resize(self%mu, nnew) + call swiftest_util_resize(self%rh, nnew) + call swiftest_util_resize(self%vh, nnew) + call swiftest_util_resize(self%rb, nnew) + call swiftest_util_resize(self%vb, nnew) + call swiftest_util_resize(self%ah, nnew) + call swiftest_util_resize(self%aobl, nnew) + call swiftest_util_resize(self%atide, nnew) + call swiftest_util_resize(self%agr, nnew) + call swiftest_util_resize(self%ir3h, nnew) + call swiftest_util_resize(self%a, nnew) + call swiftest_util_resize(self%e, nnew) + call swiftest_util_resize(self%inc, nnew) + call swiftest_util_resize(self%capom, nnew) + call swiftest_util_resize(self%omega, nnew) + call swiftest_util_resize(self%capm, nnew) self%nbody = count(self%status(1:nnew) /= INACTIVE) return @@ -2792,26 +2798,26 @@ module subroutine swiftest_util_resize_pl(self, nnew) class(swiftest_pl), intent(inout) :: self !! Swiftest massive body object integer(I4B), intent(in) :: nnew !! New size neded - call util_resize_body(self, nnew) - - call util_resize(self%mass, nnew) - call util_resize(self%Gmass, nnew) - call util_resize(self%rhill, nnew) - call util_resize(self%renc, nnew) - call util_resize(self%radius, nnew) - call util_resize(self%rbeg, nnew) - call util_resize(self%rend, nnew) - call util_resize(self%vbeg, nnew) - call util_resize(self%density, nnew) - call util_resize(self%Ip, nnew) - call util_resize(self%rot, nnew) - call util_resize(self%k2, nnew) - call util_resize(self%Q, nnew) - call util_resize(self%tlag, nnew) - call util_resize(self%kin, nnew) - call util_resize(self%lmtiny, nnew) - call util_resize(self%nplenc, nnew) - call util_resize(self%ntpenc, nnew) + call swiftest_util_resize_body(self, nnew) + + call swiftest_util_resize(self%mass, nnew) + call swiftest_util_resize(self%Gmass, nnew) + call swiftest_util_resize(self%rhill, nnew) + call swiftest_util_resize(self%renc, nnew) + call swiftest_util_resize(self%radius, nnew) + call swiftest_util_resize(self%rbeg, nnew) + call swiftest_util_resize(self%rend, nnew) + call swiftest_util_resize(self%vbeg, nnew) + call swiftest_util_resize(self%density, nnew) + call swiftest_util_resize(self%Ip, nnew) + call swiftest_util_resize(self%rot, nnew) + call swiftest_util_resize(self%k2, nnew) + call swiftest_util_resize(self%Q, nnew) + call swiftest_util_resize(self%tlag, nnew) + call swiftest_util_resize(self%kin, nnew) + call swiftest_util_resize(self%lmtiny, nnew) + call swiftest_util_resize(self%nplenc, nnew) + call swiftest_util_resize(self%ntpenc, nnew) @@ -2830,12 +2836,12 @@ module subroutine swiftest_util_resize_tp(self, nnew) class(swiftest_tp), intent(inout) :: self !! Swiftest test particle object integer(I4B), intent(in) :: nnew !! New size neded - call util_resize_body(self, nnew) + call swiftest_util_resize_body(self, nnew) - call util_resize(self%nplenc, nnew) - call util_resize(self%isperi, nnew) - call util_resize(self%peri, nnew) - call util_resize(self%atp, nnew) + call swiftest_util_resize(self%nplenc, nnew) + call swiftest_util_resize(self%isperi, nnew) + call swiftest_util_resize(self%peri, nnew) + call swiftest_util_resize(self%atp, nnew) return end subroutine swiftest_util_resize_tp @@ -3139,7 +3145,7 @@ module function swiftest_util_solve_linear_system_d(A,b,n,lerr) result(x) call ieee_set_status(original_fpe_status) return - end function util_solve_linear_system_d + end function swiftest_util_solve_linear_system_d module function swiftest_util_solve_linear_system_q(A,b,n,lerr) result(x) @@ -3175,7 +3181,7 @@ module function swiftest_util_solve_linear_system_q(A,b,n,lerr) result(x) call ieee_set_status(original_fpe_status) return - end function util_solve_linear_system_q + end function swiftest_util_solve_linear_system_q function solve_wbs(u) result(x) ! solve with backward substitution !! Based on code available on Rosetta Code: https://rosettacode.org/wiki/Gaussian_elimination#Fortran @@ -3312,7 +3318,7 @@ module function swiftest_util_solve_rkf45(f, y0in, t1, dt0, tol) result(y1) if (s >= 1.0_DP) exit ! Good step! if (i == MAXREDUX) then write(*,*) "Something has gone wrong in util_solve_rkf45!! Step size reduction has gone too far this time!" - call util_exit(FAILURE) + call swiftest_util_exit(FAILURE) end if end do @@ -3325,7 +3331,7 @@ module function swiftest_util_solve_rkf45(f, y0in, t1, dt0, tol) result(y1) end do return - end function util_solve_rkf45 + end function swiftest_util_solve_rkf45 @@ -3354,21 +3360,21 @@ module subroutine swiftest_util_sort_body(self, sortby, ascending) associate(body => self, n => self%nbody) select case(sortby) case("id") - call util_sort(direction * body%id(1:n), ind) + call swiftest_util_sort(direction * body%id(1:n), ind) case("status") - call util_sort(direction * body%status(1:n), ind) + call swiftest_util_sort(direction * body%status(1:n), ind) case("ir3h") - call util_sort(direction * body%ir3h(1:n), ind) + call swiftest_util_sort(direction * body%ir3h(1:n), ind) case("a") - call util_sort(direction * body%a(1:n), ind) + call swiftest_util_sort(direction * body%a(1:n), ind) case("e") - call util_sort(direction * body%e(1:n), ind) + call swiftest_util_sort(direction * body%e(1:n), ind) case("inc") - call util_sort(direction * body%inc(1:n), ind) + call swiftest_util_sort(direction * body%inc(1:n), ind) case("capom") - call util_sort(direction * body%capom(1:n), ind) + call swiftest_util_sort(direction * body%capom(1:n), ind) case("mu") - call util_sort(direction * body%mu(1:n), ind) + call swiftest_util_sort(direction * body%mu(1:n), ind) case("lfirst", "nbody", "ldiscard", "rh", "vh", "rb", "vb", "ah", "aobl", "atide", "agr") write(*,*) 'Cannot sort by ' // trim(adjustl(sortby)) // '. Component not sortable!' case default @@ -3393,7 +3399,7 @@ pure module subroutine swiftest_util_sort_dp(arr) ! Arguments real(DP), dimension(:), intent(inout) :: arr - call qsort_DP(arr) + call swiftest_qsort_DP(arr) return end subroutine swiftest_util_sort_dp @@ -3422,7 +3428,7 @@ pure module subroutine swiftest_util_sort_index_dp(arr, ind) end if allocate(tmparr, mold=arr) tmparr(:) = arr(ind(:)) - call qsort_DP(tmparr, ind) + call swiftest_qsort_DP(tmparr, ind) return end subroutine swiftest_util_sort_index_dp @@ -3442,13 +3448,13 @@ recursive pure subroutine swiftest_qsort_DP(arr, ind) if (size(arr) > 1) then if (present(ind)) then - call partition_DP(arr, iq, ind) - call qsort_DP(arr(:iq-1),ind(:iq-1)) - call qsort_DP(arr(iq:), ind(iq:)) + call swiftest_partition_DP(arr, iq, ind) + call swiftest_qsort_DP(arr(:iq-1),ind(:iq-1)) + call swiftest_qsort_DP(arr(iq:), ind(iq:)) else - call partition_DP(arr, iq) - call qsort_DP(arr(:iq-1)) - call qsort_DP(arr(iq:)) + call swiftest_partition_DP(arr, iq) + call swiftest_qsort_DP(arr(:iq-1)) + call swiftest_qsort_DP(arr(iq:)) end if end if @@ -3523,7 +3529,7 @@ pure module subroutine swiftest_util_sort_i4b(arr) ! Arguments integer(I4B), dimension(:), intent(inout) :: arr - call qsort_I4B(arr) + call swiftest_qsort_I4B(arr) return end subroutine swiftest_util_sort_i4b @@ -3551,7 +3557,7 @@ pure module subroutine swiftest_util_sort_index_I4B(arr, ind) end if allocate(tmparr, mold=arr) tmparr(:) = arr(ind(:)) - call qsort_I4B(tmparr, ind) + call swiftest_qsort_I4B(tmparr, ind) return end subroutine swiftest_util_sort_index_I4B @@ -3579,7 +3585,7 @@ pure module subroutine swiftest_util_sort_index_I4B_I8Bind(arr, ind) end if allocate(tmparr, mold=arr) tmparr(:) = arr(ind(:)) - call qsort_I4B_I8Bind(tmparr, ind) + call swiftest_qsort_I4B_I8Bind(tmparr, ind) return end subroutine swiftest_util_sort_index_I4B_I8Bind @@ -3599,13 +3605,13 @@ recursive pure subroutine swiftest_qsort_I4B(arr, ind) if (size(arr) > 1) then if (present(ind)) then - call partition_I4B(arr, iq, ind) - call qsort_I4B(arr(:iq-1),ind(:iq-1)) - call qsort_I4B(arr(iq:), ind(iq:)) + call swiftest_partition_I4B(arr, iq, ind) + call swiftest_qsort_I4B(arr(:iq-1),ind(:iq-1)) + call swiftest_qsort_I4B(arr(iq:), ind(iq:)) else - call partition_I4B(arr, iq) - call qsort_I4B(arr(:iq-1)) - call qsort_I4B(arr(iq:)) + call swiftest_partition_I4B(arr, iq) + call swiftest_qsort_I4B(arr(:iq-1)) + call swiftest_qsort_I4B(arr(iq:)) end if end if @@ -3626,13 +3632,13 @@ recursive pure subroutine swiftest_qsort_I4B_I8Bind(arr, ind) if (size(arr) > 1_I8B) then if (present(ind)) then - call partition_I4B_I8Bind(arr, iq, ind) - call qsort_I4B_I8Bind(arr(:iq-1_I8B),ind(:iq-1_I8B)) - call qsort_I4B_I8Bind(arr(iq:), ind(iq:)) + call swiftest_partition_I4B_I8Bind(arr, iq, ind) + call swiftest_qsort_I4B_I8Bind(arr(:iq-1_I8B),ind(:iq-1_I8B)) + call swiftest_qsort_I4B_I8Bind(arr(iq:), ind(iq:)) else - call partition_I4B_I8Bind(arr, iq) - call qsort_I4B_I8Bind(arr(:iq-1_I8B)) - call qsort_I4B_I8Bind(arr(iq:)) + call swiftest_partition_I4B_I8Bind(arr, iq) + call swiftest_qsort_I4B_I8Bind(arr(:iq-1_I8B)) + call swiftest_qsort_I4B_I8Bind(arr(iq:)) end if end if @@ -3654,13 +3660,13 @@ recursive pure subroutine swiftest_qsort_I8B_I8Bind(arr, ind) if (size(arr) > 1_I8B) then if (present(ind)) then - call partition_I8B_I8Bind(arr, iq, ind) - call qsort_I8B_I8Bind(arr(:iq-1_I8B),ind(:iq-1_I8B)) - call qsort_I8B_I8Bind(arr(iq:), ind(iq:)) + call swiftest_partition_I8B_I8Bind(arr, iq, ind) + call swiftest_qsort_I8B_I8Bind(arr(:iq-1_I8B),ind(:iq-1_I8B)) + call swiftest_qsort_I8B_I8Bind(arr(iq:), ind(iq:)) else - call partition_I8B_I8Bind(arr, iq) - call qsort_I8B_I8Bind(arr(:iq-1_I8B)) - call qsort_I8B_I8Bind(arr(iq:)) + call swiftest_partition_I8B_I8Bind(arr, iq) + call swiftest_qsort_I8B_I8Bind(arr(:iq-1_I8B)) + call swiftest_qsort_I8B_I8Bind(arr(iq:)) end if end if @@ -3846,7 +3852,7 @@ pure module subroutine swiftest_util_sort_sp(arr) ! Arguments real(SP), dimension(:), intent(inout) :: arr - call qsort_SP(arr) + call swiftest_qsort_SP(arr) return end subroutine swiftest_util_sort_sp @@ -3874,7 +3880,7 @@ pure module subroutine swiftest_util_sort_index_sp(arr, ind) end if allocate(tmparr, mold=arr) tmparr(:) = arr(ind(:)) - call qsort_SP(tmparr, ind) + call swiftest_qsort_SP(tmparr, ind) return end subroutine swiftest_util_sort_index_sp @@ -3894,13 +3900,13 @@ recursive pure subroutine swiftest_qsort_SP(arr, ind) if (size(arr) > 1) then if (present(ind)) then - call partition_SP(arr, iq, ind) - call qsort_SP(arr(:iq-1),ind(:iq-1)) - call qsort_SP(arr(iq:), ind(iq:)) + call swiftest_partition_SP(arr, iq, ind) + call swiftest_qsort_SP(arr(:iq-1),ind(:iq-1)) + call swiftest_qsort_SP(arr(iq:), ind(iq:)) else - call partition_SP(arr, iq) - call qsort_SP(arr(:iq-1)) - call qsort_SP(arr(iq:)) + call swiftest_partition_SP(arr, iq) + call swiftest_qsort_SP(arr(:iq-1)) + call swiftest_qsort_SP(arr(iq:)) end if end if @@ -3990,25 +3996,25 @@ module subroutine swiftest_util_sort_pl(self, sortby, ascending) associate(pl => self, npl => self%nbody) select case(sortby) case("Gmass","mass") - call util_sort(direction * pl%Gmass(1:npl), ind) + call swiftest_util_sort(direction * pl%Gmass(1:npl), ind) case("rhill") - call util_sort(direction * pl%rhill(1:npl), ind) + call swiftest_util_sort(direction * pl%rhill(1:npl), ind) case("renc") - call util_sort(direction * pl%renc(1:npl), ind) + call swiftest_util_sort(direction * pl%renc(1:npl), ind) case("radius") - call util_sort(direction * pl%radius(1:npl), ind) + call swiftest_util_sort(direction * pl%radius(1:npl), ind) case("density") - call util_sort(direction * pl%density(1:npl), ind) + call swiftest_util_sort(direction * pl%density(1:npl), ind) case("k2") - call util_sort(direction * pl%k2(1:npl), ind) + call swiftest_util_sort(direction * pl%k2(1:npl), ind) case("Q") - call util_sort(direction * pl%Q(1:npl), ind) + call swiftest_util_sort(direction * pl%Q(1:npl), ind) case("tlag") - call util_sort(direction * pl%tlag(1:npl), ind) + call swiftest_util_sort(direction * pl%tlag(1:npl), ind) case("rbeg", "rend", "vbeg", "Ip", "rot", "k_plpl", "nplpl") write(*,*) 'Cannot sort by ' // trim(adjustl(sortby)) // '. Component not sortable!' case default ! Look for components in the parent class - call util_sort_body(pl, sortby, ascending) + call swiftest_util_sort_body(pl, sortby, ascending) return end select @@ -4045,13 +4051,13 @@ module subroutine swiftest_util_sort_tp(self, sortby, ascending) associate(tp => self, ntp => self%nbody) select case(sortby) case("peri") - call util_sort(direction * tp%peri(1:ntp), ind) + call swiftest_util_sort(direction * tp%peri(1:ntp), ind) case("atp") - call util_sort(direction * tp%atp(1:ntp), ind) + call swiftest_util_sort(direction * tp%atp(1:ntp), ind) case("isperi") write(*,*) 'Cannot sort by ' // trim(adjustl(sortby)) // '. Component not sortable!' case default ! Look for components in the parent class - call util_sort_body(tp, sortby, ascending) + call swiftest_util_sort_body(tp, sortby, ascending) return end select @@ -4074,32 +4080,32 @@ module subroutine swiftest_util_sort_rearrange_body(self, ind) integer(I4B), dimension(:), intent(in) :: ind !! Index array used to restructure the body (should contain all 1:n index values in the desired order) associate(n => self%nbody) - call util_sort_rearrange(self%id, ind, n) - call util_sort_rearrange(self%lmask, ind, n) - call util_sort_rearrange(self%info, ind, n) - call util_sort_rearrange(self%status, ind, n) - call util_sort_rearrange(self%ldiscard, ind, n) - call util_sort_rearrange(pl%lcollision, ind, n) - call util_sort_rearrange(pl%lencounter, ind, n) - call util_sort_rearrange(self%rh, ind, n) - call util_sort_rearrange(self%vh, ind, n) - call util_sort_rearrange(self%rb, ind, n) - call util_sort_rearrange(self%vb, ind, n) - call util_sort_rearrange(self%ah, ind, n) - call util_sort_rearrange(self%aobl, ind, n) - call util_sort_rearrange(self%agr, ind, n) - call util_sort_rearrange(self%atide, ind, n) - call util_sort_rearrange(self%ir3h, ind, n) - call util_sort_rearrange(self%isperi, ind, n) - call util_sort_rearrange(self%peri, ind, n) - call util_sort_rearrange(self%atp, ind, n) - call util_sort_rearrange(self%mu, ind, n) - call util_sort_rearrange(self%a, ind, n) - call util_sort_rearrange(self%e, ind, n) - call util_sort_rearrange(self%inc, ind, n) - call util_sort_rearrange(self%capom, ind, n) - call util_sort_rearrange(self%omega, ind, n) - call util_sort_rearrange(self%capm, ind, n) + call swiftest_util_sort_rearrange(self%id, ind, n) + call swiftest_util_sort_rearrange(self%lmask, ind, n) + call swiftest_util_sort_rearrange(self%info, ind, n) + call swiftest_util_sort_rearrange(self%status, ind, n) + call swiftest_util_sort_rearrange(self%ldiscard, ind, n) + call swiftest_util_sort_rearrange(self%lcollision, ind, n) + call swiftest_util_sort_rearrange(self%lencounter, ind, n) + call swiftest_util_sort_rearrange(self%rh, ind, n) + call swiftest_util_sort_rearrange(self%vh, ind, n) + call swiftest_util_sort_rearrange(self%rb, ind, n) + call swiftest_util_sort_rearrange(self%vb, ind, n) + call swiftest_util_sort_rearrange(self%ah, ind, n) + call swiftest_util_sort_rearrange(self%aobl, ind, n) + call swiftest_util_sort_rearrange(self%agr, ind, n) + call swiftest_util_sort_rearrange(self%atide, ind, n) + call swiftest_util_sort_rearrange(self%ir3h, ind, n) + call swiftest_util_sort_rearrange(self%isperi, ind, n) + call swiftest_util_sort_rearrange(self%peri, ind, n) + call swiftest_util_sort_rearrange(self%atp, ind, n) + call swiftest_util_sort_rearrange(self%mu, ind, n) + call swiftest_util_sort_rearrange(self%a, ind, n) + call swiftest_util_sort_rearrange(self%e, ind, n) + call swiftest_util_sort_rearrange(self%inc, ind, n) + call swiftest_util_sort_rearrange(self%capom, ind, n) + call swiftest_util_sort_rearrange(self%omega, ind, n) + call swiftest_util_sort_rearrange(self%capm, ind, n) end associate return @@ -4267,7 +4273,7 @@ module subroutine swiftest_util_sort_rearrange_arr_info(arr, ind, n) if (.not. allocated(arr) .or. n <= 0) return allocate(tmp, mold=arr) - call util_copy_particle_info_arr(arr, tmp, ind) + call swiftest_util_copy_particle_info_arr(arr, tmp, ind) call move_alloc(tmp, arr) return @@ -4284,27 +4290,27 @@ module subroutine swiftest_util_sort_rearrange_pl(self, ind) integer(I4B), dimension(:), intent(in) :: ind !! Index array used to restructure the body (should contain all 1:n index values in the desired order) associate(pl => self, npl => self%nbody) - call util_sort_rearrange(pl%mass, ind, npl) - call util_sort_rearrange(pl%Gmass, ind, npl) - call util_sort_rearrange(pl%rhill, ind, npl) - call util_sort_rearrange(pl%renc, ind, npl) - call util_sort_rearrange(pl%radius, ind, npl) - call util_sort_rearrange(pl%density, ind, npl) - call util_sort_rearrange(pl%rbeg, ind, npl) - call util_sort_rearrange(pl%vbeg, ind, npl) - call util_sort_rearrange(pl%Ip, ind, npl) - call util_sort_rearrange(pl%rot, ind, npl) - call util_sort_rearrange(pl%k2, ind, npl) - call util_sort_rearrange(pl%Q, ind, npl) - call util_sort_rearrange(pl%tlag, ind, npl) - call util_sort_rearrange(pl%kin, ind, npl) - call util_sort_rearrange(pl%lmtiny, ind, npl) - call util_sort_rearrange(pl%nplenc, ind, npl) - call util_sort_rearrange(pl%ntpenc, ind, npl) + call swiftest_util_sort_rearrange(pl%mass, ind, npl) + call swiftest_util_sort_rearrange(pl%Gmass, ind, npl) + call swiftest_util_sort_rearrange(pl%rhill, ind, npl) + call swiftest_util_sort_rearrange(pl%renc, ind, npl) + call swiftest_util_sort_rearrange(pl%radius, ind, npl) + call swiftest_util_sort_rearrange(pl%density, ind, npl) + call swiftest_util_sort_rearrange(pl%rbeg, ind, npl) + call swiftest_util_sort_rearrange(pl%vbeg, ind, npl) + call swiftest_util_sort_rearrange(pl%Ip, ind, npl) + call swiftest_util_sort_rearrange(pl%rot, ind, npl) + call swiftest_util_sort_rearrange(pl%k2, ind, npl) + call swiftest_util_sort_rearrange(pl%Q, ind, npl) + call swiftest_util_sort_rearrange(pl%tlag, ind, npl) + call swiftest_util_sort_rearrange(pl%kin, ind, npl) + call swiftest_util_sort_rearrange(pl%lmtiny, ind, npl) + call swiftest_util_sort_rearrange(pl%nplenc, ind, npl) + call swiftest_util_sort_rearrange(pl%ntpenc, ind, npl) if (allocated(pl%k_plpl)) deallocate(pl%k_plpl) - call util_sort_rearrange_body(pl, ind) + call swiftest_util_sort_rearrange_body(pl, ind) end associate return @@ -4322,11 +4328,11 @@ module subroutine swiftest_util_sort_rearrange_tp(self, ind) integer(I4B), dimension(:), intent(in) :: ind !! Index array used to restructure the body (should contain all 1:n index values in the desired order) associate(tp => self, ntp => self%nbody) - call util_sort_rearrange(tp%nplenc, ind, ntp) + call swiftest_util_sort_rearrange(tp%nplenc, ind, ntp) if (allocated(tp%k_pltp)) deallocate(tp%k_pltp) - call util_sort_rearrange_body(tp, ind) + call swiftest_util_sort_rearrange_body(tp, ind) end associate return @@ -4577,14 +4583,14 @@ module subroutine swiftest_util_spill_arr_info(keeps, discards, lspill_list, lde allocate(idx(nspill)) idx(:) = pack([(i, i = 1, nlist)], lspill_list) - call util_copy_particle_info_arr(keeps, discards, idx) + call swiftest_util_copy_particle_info_arr(keeps, discards, idx) if (ldestructive) then if (nkeep > 0) then deallocate(idx) allocate(idx(nkeep)) allocate(tmp(nkeep)) idx(:) = pack([(i, i = 1, nlist)], .not. lspill_list) - call util_copy_particle_info_arr(keeps, tmp, idx) + call swiftest_util_copy_particle_info_arr(keeps, tmp, idx) call move_alloc(tmp, keeps) else deallocate(keeps) @@ -4655,32 +4661,32 @@ module subroutine swiftest_util_spill_body(self, discards, lspill_list, ldestruc !> Spill all the common components associate(keeps => self) - call util_spill(keeps%id, discards%id, lspill_list, ldestructive) - call util_spill(keeps%info, discards%info, lspill_list, ldestructive) - call util_spill(keeps%lmask, discards%lmask, lspill_list, ldestructive) - call util_spill(keeps%status, discards%status, lspill_list, ldestructive) - call util_spill(keeps%ldiscard, discards%ldiscard, lspill_list, ldestructive) - call util_spill(keeps%lcollision, discards%lcollision, lspill_list, ldestructive) - call util_spill(keeps%lencounter, discards%lencounter, lspill_list, ldestructive) - call util_spill(keeps%mu, discards%mu, lspill_list, ldestructive) - call util_spill(keeps%rh, discards%rh, lspill_list, ldestructive) - call util_spill(keeps%vh, discards%vh, lspill_list, ldestructive) - call util_spill(keeps%rb, discards%rb, lspill_list, ldestructive) - call util_spill(keeps%vb, discards%vb, lspill_list, ldestructive) - call util_spill(keeps%ah, discards%ah, lspill_list, ldestructive) - call util_spill(keeps%aobl, discards%aobl, lspill_list, ldestructive) - call util_spill(keeps%agr, discards%agr, lspill_list, ldestructive) - call util_spill(keeps%atide, discards%atide, lspill_list, ldestructive) - call util_spill(keeps%ir3h, discards%ir3h, lspill_list, ldestructive) - call util_spill(keeps%isperi, discards%isperi, lspill_list, ldestructive) - call util_spill(keeps%peri, discards%peri, lspill_list, ldestructive) - call util_spill(keeps%atp, discards%atp, lspill_list, ldestructive) - call util_spill(keeps%a, discards%a, lspill_list, ldestructive) - call util_spill(keeps%e, discards%e, lspill_list, ldestructive) - call util_spill(keeps%inc, discards%inc, lspill_list, ldestructive) - call util_spill(keeps%capom, discards%capom, lspill_list, ldestructive) - call util_spill(keeps%omega, discards%omega, lspill_list, ldestructive) - call util_spill(keeps%capm, discards%capm, lspill_list, ldestructive) + call swiftest_util_spill(keeps%id, discards%id, lspill_list, ldestructive) + call swiftest_util_spill(keeps%info, discards%info, lspill_list, ldestructive) + call swiftest_util_spill(keeps%lmask, discards%lmask, lspill_list, ldestructive) + call swiftest_util_spill(keeps%status, discards%status, lspill_list, ldestructive) + call swiftest_util_spill(keeps%ldiscard, discards%ldiscard, lspill_list, ldestructive) + call swiftest_util_spill(keeps%lcollision, discards%lcollision, lspill_list, ldestructive) + call swiftest_util_spill(keeps%lencounter, discards%lencounter, lspill_list, ldestructive) + call swiftest_util_spill(keeps%mu, discards%mu, lspill_list, ldestructive) + call swiftest_util_spill(keeps%rh, discards%rh, lspill_list, ldestructive) + call swiftest_util_spill(keeps%vh, discards%vh, lspill_list, ldestructive) + call swiftest_util_spill(keeps%rb, discards%rb, lspill_list, ldestructive) + call swiftest_util_spill(keeps%vb, discards%vb, lspill_list, ldestructive) + call swiftest_util_spill(keeps%ah, discards%ah, lspill_list, ldestructive) + call swiftest_util_spill(keeps%aobl, discards%aobl, lspill_list, ldestructive) + call swiftest_util_spill(keeps%agr, discards%agr, lspill_list, ldestructive) + call swiftest_util_spill(keeps%atide, discards%atide, lspill_list, ldestructive) + call swiftest_util_spill(keeps%ir3h, discards%ir3h, lspill_list, ldestructive) + call swiftest_util_spill(keeps%isperi, discards%isperi, lspill_list, ldestructive) + call swiftest_util_spill(keeps%peri, discards%peri, lspill_list, ldestructive) + call swiftest_util_spill(keeps%atp, discards%atp, lspill_list, ldestructive) + call swiftest_util_spill(keeps%a, discards%a, lspill_list, ldestructive) + call swiftest_util_spill(keeps%e, discards%e, lspill_list, ldestructive) + call swiftest_util_spill(keeps%inc, discards%inc, lspill_list, ldestructive) + call swiftest_util_spill(keeps%capom, discards%capom, lspill_list, ldestructive) + call swiftest_util_spill(keeps%omega, discards%omega, lspill_list, ldestructive) + call swiftest_util_spill(keeps%capm, discards%capm, lspill_list, ldestructive) nbody_old = keeps%nbody @@ -4710,28 +4716,28 @@ module subroutine swiftest_util_spill_pl(self, discards, lspill_list, ldestructi select type (discards) ! The standard requires us to select the type of both arguments in order to access all the components class is (swiftest_pl) !> Spill components specific to the massive body class - call util_spill(keeps%mass, discards%mass, lspill_list, ldestructive) - call util_spill(keeps%Gmass, discards%Gmass, lspill_list, ldestructive) - call util_spill(keeps%rhill, discards%rhill, lspill_list, ldestructive) - call util_spill(keeps%renc, discards%renc, lspill_list, ldestructive) - call util_spill(keeps%radius, discards%radius, lspill_list, ldestructive) - call util_spill(keeps%density, discards%density, lspill_list, ldestructive) - call util_spill(keeps%rbeg, discards%rbeg, lspill_list, ldestructive) - call util_spill(keeps%rend, discards%rend, lspill_list, ldestructive) - call util_spill(keeps%vbeg, discards%vbeg, lspill_list, ldestructive) - call util_spill(keeps%Ip, discards%Ip, lspill_list, ldestructive) - call util_spill(keeps%rot, discards%rot, lspill_list, ldestructive) - call util_spill(keeps%k2, discards%k2, lspill_list, ldestructive) - call util_spill(keeps%Q, discards%Q, lspill_list, ldestructive) - call util_spill(keeps%tlag, discards%tlag, lspill_list, ldestructive) - call util_spill(keeps%kin, discards%kin, lspill_list, ldestructive) - call util_spill(keeps%lmtiny, discards%lmtiny, lspill_list, ldestructive) - call util_spill(keeps%nplenc, discards%nplenc, lspill_list, ldestructive) - call util_spill(keeps%ntpenc, discards%ntpenc, lspill_list, ldestructive) + call swiftest_util_spill(keeps%mass, discards%mass, lspill_list, ldestructive) + call swiftest_util_spill(keeps%Gmass, discards%Gmass, lspill_list, ldestructive) + call swiftest_util_spill(keeps%rhill, discards%rhill, lspill_list, ldestructive) + call swiftest_util_spill(keeps%renc, discards%renc, lspill_list, ldestructive) + call swiftest_util_spill(keeps%radius, discards%radius, lspill_list, ldestructive) + call swiftest_util_spill(keeps%density, discards%density, lspill_list, ldestructive) + call swiftest_util_spill(keeps%rbeg, discards%rbeg, lspill_list, ldestructive) + call swiftest_util_spill(keeps%rend, discards%rend, lspill_list, ldestructive) + call swiftest_util_spill(keeps%vbeg, discards%vbeg, lspill_list, ldestructive) + call swiftest_util_spill(keeps%Ip, discards%Ip, lspill_list, ldestructive) + call swiftest_util_spill(keeps%rot, discards%rot, lspill_list, ldestructive) + call swiftest_util_spill(keeps%k2, discards%k2, lspill_list, ldestructive) + call swiftest_util_spill(keeps%Q, discards%Q, lspill_list, ldestructive) + call swiftest_util_spill(keeps%tlag, discards%tlag, lspill_list, ldestructive) + call swiftest_util_spill(keeps%kin, discards%kin, lspill_list, ldestructive) + call swiftest_util_spill(keeps%lmtiny, discards%lmtiny, lspill_list, ldestructive) + call swiftest_util_spill(keeps%nplenc, discards%nplenc, lspill_list, ldestructive) + call swiftest_util_spill(keeps%ntpenc, discards%ntpenc, lspill_list, ldestructive) if (ldestructive .and. allocated(keeps%k_plpl)) deallocate(keeps%k_plpl) - call util_spill_body(keeps, discards, lspill_list, ldestructive) + call swiftest_util_spill_body(keeps, discards, lspill_list, ldestructive) class default write(*,*) 'Error! spill method called for incompatible return type on swiftest_pl' end select @@ -4757,8 +4763,8 @@ module subroutine swiftest_util_spill_tp(self, discards, lspill_list, ldestructi select type(discards) class is (swiftest_tp) !> Spill components specific to the test particle class - call util_spill(keeps%nplenc, discards%nplenc, lspill_list, ldestructive) - call util_spill_body(keeps, discards, lspill_list, ldestructive) + call swiftest_util_spill(keeps%nplenc, discards%nplenc, lspill_list, ldestructive) + call swiftest_util_spill_body(keeps, discards, lspill_list, ldestructive) class default write(*,*) 'Error! spill method called for incompatible return type on swiftest_tp' end select @@ -4857,12 +4863,12 @@ module subroutine swiftest_util_valid_id_system(self, param) do i = 1, ntp idarr(1+npl+i) = tp%id(i) end do - call util_sort(idarr) + call swiftest_util_sort(idarr) do i = 1, npl + ntp if (idarr(i) == idarr(i+1)) then write(*, *) "Swiftest error:" write(*, *) " more than one body/particle has id = ", idarr(i) - call util_exit(FAILURE) + call swiftest_util_exit(FAILURE) end if end do param%maxid = max(param%maxid, maxval(idarr)) diff --git a/src/symba/symba_discard.f90 b/src/symba/symba_discard.f90 index a14911daa..5cdf5eab5 100644 --- a/src/symba/symba_discard.f90 +++ b/src/symba/symba_discard.f90 @@ -225,26 +225,23 @@ subroutine symba_discard_nonplpl(pl, system, param) class(symba_pl), allocatable :: plsub ! First check for collisions with the central body - associate(npl => pl%nbody, cb => system%cb) + associate(npl => pl%nbody, cb => system%cb, pl_discards => system%pl_discards) if (npl == 0) return - select type(pl_discards => system%pl_discards) - class is (symba_merger) - if ((param%rmin >= 0.0_DP) .or. (param%rmax >= 0.0_DP) .or. (param%rmaxu >= 0.0_DP)) then - call symba_discard_cb_pl(pl, system, param) - end if - if (param%qmin >= 0.0_DP) call symba_discard_peri_pl(pl, system, param) - if (any(pl%ldiscard(1:npl))) then - ldiscard(1:npl) = pl%ldiscard(1:npl) - - allocate(plsub, mold=pl) - call pl%spill(plsub, ldiscard, ldestructive=.false.) - nsub = plsub%nbody - nstart = pl_discards%nbody + 1 - nend = pl_discards%nbody + nsub - call pl_discards%append(plsub, lsource_mask=[(.true., i = 1, nsub)]) - - end if - end select + if ((param%rmin >= 0.0_DP) .or. (param%rmax >= 0.0_DP) .or. (param%rmaxu >= 0.0_DP)) then + call symba_discard_cb_pl(pl, system, param) + end if + if (param%qmin >= 0.0_DP) call symba_discard_peri_pl(pl, system, param) + if (any(pl%ldiscard(1:npl))) then + ldiscard(1:npl) = pl%ldiscard(1:npl) + + allocate(plsub, mold=pl) + call pl%spill(plsub, ldiscard, ldestructive=.false.) + nsub = plsub%nbody + nstart = pl_discards%nbody + 1 + nend = pl_discards%nbody + nsub + call pl_discards%append(plsub, lsource_mask=[(.true., i = 1, nsub)]) + + end if end associate return diff --git a/src/symba/symba_encounter_check.f90 b/src/symba/symba_encounter_check.f90 index eb05f0664..d0fc96990 100644 --- a/src/symba/symba_encounter_check.f90 +++ b/src/symba/symba_encounter_check.f90 @@ -87,16 +87,9 @@ module function symba_encounter_check_pl(self, param, system, dt, irec) result(l end function symba_encounter_check_pl - module function symba_encounter_io_netcdf_check(self, param, system, dt, irec) result(lany_encounter) - !! author: David A. Minton - !! - !! Check for an encounter between test particles and massive bodies in the plplenc and pltpenc list. - !! Note: This method works for the polymorphic pltp_encounter and plpl_encounter types. - !! - !! Adapted from portions of David E. Kaufmann's Swifter routine: symba_step_recur.f90 + module function symba_encounter_check_list_plpl(self, param, system, dt, irec) result(lany_encounter) implicit none - ! Arguments - class(symba_encounter), intent(inout) :: self !! SyMBA pl-pl encounter list object + class(symba_list_plpl), intent(inout) :: self !! SyMBA pl-pl encounter list object class(swiftest_parameters), intent(inout) :: param !! Current swiftest run configuration parameters class(symba_nbody_system), intent(inout) :: system !! SyMBA nbody system object real(DP), intent(in) :: dt !! step size @@ -105,7 +98,6 @@ module function symba_encounter_io_netcdf_check(self, param, system, dt, irec) r ! Internals integer(I4B) :: i, j, k, lidx, nenc_enc real(DP), dimension(NDIM) :: xr, vr - logical :: isplpl real(DP) :: rlim2, rji2, rcrit12 logical, dimension(:), allocatable :: lencmask, lencounter integer(I4B), dimension(:), allocatable :: eidx @@ -113,84 +105,125 @@ module function symba_encounter_io_netcdf_check(self, param, system, dt, irec) r lany_encounter = .false. if (self%nenc == 0) return - select type(self) - class is (plpl_encounter) - isplpl = .true. - class is (pltp_encounter) - isplpl = .false. + select type(pl => system%pl) + class is (symba_pl) + allocate(lencmask(self%nenc)) + lencmask(:) = (self%status(1:self%nenc) == ACTIVE) .and. (self%level(1:self%nenc) == irec - 1) + nenc_enc = count(lencmask(:)) + if (nenc_enc == 0) return + + call pl%set_renc(irec) + + allocate(eidx(nenc_enc)) + allocate(lencounter(nenc_enc)) + eidx(:) = pack([(k, k = 1, self%nenc)], lencmask(:)) + lencounter(:) = .false. + + do concurrent(lidx = 1:nenc_enc) + k = eidx(lidx) + i = self%index1(k) + j = self%index2(k) + xr(:) = pl%rh(:,j) - pl%rh(:,i) + vr(:) = pl%vb(:,j) - pl%vb(:,i) + rcrit12 = pl%renc(i) + pl%renc(j) + call encounter_check_one(xr(1), xr(2), xr(3), vr(1), vr(2), vr(3), rcrit12, dt, lencounter(lidx), self%lvdotr(k)) + if (lencounter(lidx)) then + rlim2 = (pl%radius(i) + pl%radius(j))**2 + rji2 = dot_product(xr(:), xr(:))! Check to see if these are physically overlapping bodies first, which we should ignore + lencounter(lidx) = rji2 > rlim2 + end if + end do + + lany_encounter = any(lencounter(:)) + if (lany_encounter) then + nenc_enc = count(lencounter(:)) + eidx(1:nenc_enc) = pack(eidx(:), lencounter(:)) + do lidx = 1, nenc_enc + k = eidx(lidx) + i = self%index1(k) + j = self%index2(k) + pl%levelg(i) = irec + pl%levelm(i) = MAX(irec, pl%levelm(i)) + pl%levelg(j) = irec + pl%levelm(j) = MAX(irec, pl%levelm(j)) + self%level(k) = irec + end do + end if end select + return + end function symba_encounter_check_list_plpl + + + module function symba_encounter_check_list_pltp(self, param, system, dt, irec) result(lany_encounter) + implicit none + class(symba_list_pltp), intent(inout) :: self !! SyMBA pl-tp encounter list object + class(swiftest_parameters), intent(inout) :: param !! Current swiftest run configuration parameters + class(symba_nbody_system), intent(inout) :: system !! SyMBA nbody system object + real(DP), intent(in) :: dt !! step size + integer(I4B), intent(in) :: irec !! Current recursion level + logical :: lany_encounter !! Returns true if there is at least one close encounter + ! Internals + integer(I4B) :: i, j, k, lidx, nenc_enc + real(DP), dimension(NDIM) :: xr, vr + real(DP) :: rlim2, rji2 + logical, dimension(:), allocatable :: lencmask, lencounter + integer(I4B), dimension(:), allocatable :: eidx + + lany_encounter = .false. + if (self%nenc == 0) return + select type(pl => system%pl) class is (symba_pl) - select type(tp => system%tp) - class is (symba_tp) - allocate(lencmask(self%nenc)) - lencmask(:) = (self%status(1:self%nenc) == ACTIVE) .and. (self%level(1:self%nenc) == irec - 1) - nenc_enc = count(lencmask(:)) - if (nenc_enc == 0) return - - call pl%set_renc(irec) - - allocate(eidx(nenc_enc)) - allocate(lencounter(nenc_enc)) - eidx(:) = pack([(k, k = 1, self%nenc)], lencmask(:)) - lencounter(:) = .false. - if (isplpl) then - do concurrent(lidx = 1:nenc_enc) - k = eidx(lidx) - i = self%index1(k) - j = self%index2(k) - xr(:) = pl%rh(:,j) - pl%rh(:,i) - vr(:) = pl%vb(:,j) - pl%vb(:,i) - rcrit12 = pl%renc(i) + pl%renc(j) - call encounter_check_one(xr(1), xr(2), xr(3), vr(1), vr(2), vr(3), rcrit12, dt, lencounter(lidx), self%lvdotr(k)) - if (lencounter(lidx)) then - rlim2 = (pl%radius(i) + pl%radius(j))**2 - rji2 = dot_product(xr(:), xr(:))! Check to see if these are physically overlapping bodies first, which we should ignore - lencounter(lidx) = rji2 > rlim2 - end if - end do - else - do concurrent(lidx = 1:nenc_enc) - k = eidx(lidx) - i = self%index1(k) - j = self%index2(k) - xr(:) = tp%rh(:,j) - pl%rh(:,i) - vr(:) = tp%vb(:,j) - pl%vb(:,i) - call encounter_check_one(xr(1), xr(2), xr(3), vr(1), vr(2), vr(3), pl%renc(i), dt, & - lencounter(lidx), self%lvdotr(k)) - if (lencounter(lidx)) then - rlim2 = (pl%radius(i))**2 - rji2 = dot_product(xr(:), xr(:))! Check to see if these are physically overlapping bodies first, which we should ignore - lencounter(lidx) = rji2 > rlim2 - end if - end do + select type(tp => system%tp) + class is (symba_tp) + allocate(lencmask(self%nenc)) + lencmask(:) = (self%status(1:self%nenc) == ACTIVE) .and. (self%level(1:self%nenc) == irec - 1) + nenc_enc = count(lencmask(:)) + if (nenc_enc == 0) return + + call pl%set_renc(irec) + + allocate(eidx(nenc_enc)) + allocate(lencounter(nenc_enc)) + eidx(:) = pack([(k, k = 1, self%nenc)], lencmask(:)) + lencounter(:) = .false. + + do concurrent(lidx = 1:nenc_enc) + k = eidx(lidx) + i = self%index1(k) + j = self%index2(k) + xr(:) = tp%rh(:,j) - pl%rh(:,i) + vr(:) = tp%vb(:,j) - pl%vb(:,i) + call encounter_check_one(xr(1), xr(2), xr(3), vr(1), vr(2), vr(3), pl%renc(i), dt, & + lencounter(lidx), self%lvdotr(k)) + if (lencounter(lidx)) then + rlim2 = (pl%radius(i))**2 + rji2 = dot_product(xr(:), xr(:))! Check to see if these are physically overlapping bodies first, which we should ignore + lencounter(lidx) = rji2 > rlim2 end if - lany_encounter = any(lencounter(:)) - if (lany_encounter) then - nenc_enc = count(lencounter(:)) - eidx(1:nenc_enc) = pack(eidx(:), lencounter(:)) - do lidx = 1, nenc_enc - k = eidx(lidx) - i = self%index1(k) - j = self%index2(k) - pl%levelg(i) = irec - pl%levelm(i) = MAX(irec, pl%levelm(i)) - if (isplpl) then - pl%levelg(j) = irec - pl%levelm(j) = MAX(irec, pl%levelm(j)) - else - tp%levelg(j) = irec - tp%levelm(j) = MAX(irec, tp%levelm(j)) - end if - self%level(k) = irec - end do - end if - end select + end do + + lany_encounter = any(lencounter(:)) + if (lany_encounter) then + nenc_enc = count(lencounter(:)) + eidx(1:nenc_enc) = pack(eidx(:), lencounter(:)) + do lidx = 1, nenc_enc + k = eidx(lidx) + i = self%index1(k) + j = self%index2(k) + pl%levelg(i) = irec + pl%levelm(i) = MAX(irec, pl%levelm(i)) + tp%levelg(j) = irec + tp%levelm(j) = MAX(irec, tp%levelm(j)) + self%level(k) = irec + end do + end if + end select end select return - end function symba_encounter_check + end function symba_encounter_check_list_pltp module function symba_encounter_check_tp(self, param, system, dt, irec) result(lany_encounter) diff --git a/src/symba/symba_kick.f90 b/src/symba/symba_kick.f90 index 03b72c143..80d3ccbc2 100644 --- a/src/symba/symba_kick.f90 +++ b/src/symba/symba_kick.f90 @@ -23,33 +23,33 @@ module subroutine symba_kick_getacch_int_pl(self, param) class(symba_pl), intent(inout) :: self !! SyMBA massive body object class(swiftest_parameters), intent(inout) :: param !! Current swiftest run configuration parameter ! Internals - type(interaction_timer), save :: itimer - logical, save :: lfirst = .true. - - if (param%ladaptive_interactions) then - if (self%nplplm > 0) then - if (lfirst) then - write(itimer%loopname, *) "symba_kick_getacch_int_pl" - write(itimer%looptype, *) "INTERACTION" - call itimer%time_this_loop(param, self%nplplm, self) - lfirst = .false. - else - if (itimer%io_netcdf_check(param, self%nplplm)) call itimer%time_this_loop(param, self%nplplm, self) - end if - else - param%lflatten_interactions = .false. - end if - end if - - if (param%lflatten_interactions) then - call swiftest_kick_getacch_int_all_flat_pl(self%nbody, self%nplplm, self%k_plpl, self%rh, self%Gmass, self%radius, self%ah) - else + ! type(interaction_timer), save :: itimer + ! logical, save :: lfirst = .true. + + ! if (param%ladaptive_interactions) then + ! if (self%nplplm > 0) then + ! if (lfirst) then + ! write(itimer%loopname, *) "symba_kick_getacch_int_pl" + ! write(itimer%looptype, *) "INTERACTION" + ! call itimer%time_this_loop(param, self%nplplm, self) + ! lfirst = .false. + ! else + ! if (itimer%io_netcdf_check(param, self%nplplm)) call itimer%time_this_loop(param, self%nplplm, self) + ! end if + ! else + ! param%lflatten_interactions = .false. + ! end if + ! end if + + ! if (param%lflatten_interactions) then + ! call swiftest_kick_getacch_int_all_flat_pl(self%nbody, self%nplplm, self%k_plpl, self%rh, self%Gmass, self%radius, self%ah) + ! else call swiftest_kick_getacch_int_all_triangular_pl(self%nbody, self%nplm, self%rh, self%Gmass, self%radius, self%ah) - end if + ! end if - if (param%ladaptive_interactions .and. self%nplplm > 0) then - if (itimer%is_on) call itimer%adapt(param, self%nplplm, self) - end if + ! if (param%ladaptive_interactions .and. self%nplplm > 0) then + ! if (itimer%is_on) call itimer%adapt(param, self%nplplm, self) + ! end if return end subroutine symba_kick_getacch_int_pl @@ -144,17 +144,16 @@ module subroutine symba_kick_getacch_tp(self, system, param, t, lbeg) end subroutine symba_kick_getacch_tp - module subroutine symba_kick_encounter(self, system, dt, irec, sgn) + module subroutine symba_kick_list_plpl(self, system, dt, irec, sgn) !! author: David A. Minton !! - !! Kick barycentric velocities of massive bodies and ACTIVE test particles within SyMBA recursion. - !! Note: This method works for the polymorphic pltp_encounter and plpl_encounter types + !! Kick barycentric velocities of massive bodies within SyMBA recursion. !! !! Adapted from David E. Kaufmann's Swifter routine: symba_kick.f90 !! Adapted from Hal Levison's Swift routine symba5_kick.f implicit none ! Arguments - class(symba_encounter), intent(in) :: self !! SyMBA pl-tp encounter list object + class(symba_list_plpl), intent(in) :: self !! SyMBA pl-tp encounter list object class(symba_nbody_system), intent(inout) :: system !! SyMBA nbody system object real(DP), intent(in) :: dt !! step size integer(I4B), intent(in) :: irec !! Current recursion level @@ -164,130 +163,195 @@ module subroutine symba_kick_encounter(self, system, dt, irec, sgn) integer(I8B) :: k real(DP) :: r, rr, ri, ris, rim1, r2, ir3, fac, faci, facj real(DP), dimension(NDIM) :: dx - logical :: isplpl logical, dimension(:), allocatable :: lgoodlevel integer(I4B), dimension(:), allocatable :: good_idx if (self%nenc == 0) return - select type(self) - class is (plpl_encounter) - isplpl = .true. - class is (pltp_encounter) - isplpl = .false. - end select select type(pl => system%pl) class is (symba_pl) - select type(tp => system%tp) - class is (symba_tp) - associate(npl => pl%nbody, ntp => tp%nbody, nenc => self%nenc) - if (npl == 0) return - pl%lmask(1:npl) = pl%status(1:npl) /= INACTIVE - if (.not. isplpl) then - if (ntp == 0) return - tp%lmask(1:ntp) = tp%status(1:ntp) /= INACTIVE - end if - allocate(lgoodlevel(nenc)) + associate(npl => pl%nbody, nenc => self%nenc) + if (npl == 0) return + pl%lmask(1:npl) = pl%status(1:npl) /= INACTIVE + allocate(lgoodlevel(nenc)) - irm1 = irec - 1 + irm1 = irec - 1 - if (sgn < 0) then - irecl = irec - 1 - else - irecl = irec - end if + if (sgn < 0) then + irecl = irec - 1 + else + irecl = irec + end if - do k = 1, nenc - i = self%index1(k) - j = self%index2(k) - if (isplpl) then - lgoodlevel(k) = (pl%levelg(i) >= irm1) .and. (pl%levelg(j) >= irm1) + do k = 1, nenc + i = self%index1(k) + j = self%index2(k) + lgoodlevel(k) = (pl%levelg(i) >= irm1) .and. (pl%levelg(j) >= irm1) + lgoodlevel(k) = (self%status(k) == ACTIVE) .and. lgoodlevel(k) + end do + ngood = count(lgoodlevel(:)) + if (ngood > 0_I8B) then + allocate(good_idx(ngood)) + good_idx(:) = pack([(i, i = 1, nenc)], lgoodlevel(:)) + + do concurrent (k = 1:ngood) + i = self%index1(good_idx(k)) + j = self%index2(good_idx(k)) + pl%ah(:,i) = 0.0_DP + pl%ah(:,j) = 0.0_DP + end do + + do k = 1, ngood + i = self%index1(good_idx(k)) + j = self%index2(good_idx(k)) + ri = ((pl%rhill(i) + pl%rhill(j))**2) * (RHSCALE**2) * (RSHELL**(2*irecl)) + rim1 = ri * (RSHELL**2) + dx(:) = pl%rh(:,j) - pl%rh(:,i) + + r2 = dot_product(dx(:), dx(:)) + if (r2 < rim1) then + fac = 0.0_DP + lgoodlevel(good_idx(k)) = .false. + cycle + end if + if (r2 < ri) then + ris = sqrt(ri) + r = sqrt(r2) + rr = (ris - r) / (ris * (1.0_DP - RSHELL)) + fac = (r2**(-1.5_DP)) * (1.0_DP - 3 * (rr**2) + 2 * (rr**3)) else - lgoodlevel(k) = (pl%levelg(i) >= irm1) .and. (tp%levelg(j) >= irm1) + ir3 = 1.0_DP / (r2 * sqrt(r2)) + fac = ir3 end if - lgoodlevel(k) = (self%status(k) == ACTIVE) .and. lgoodlevel(k) + faci = fac * pl%Gmass(i) + facj = fac * pl%Gmass(j) + pl%ah(:, i) = pl%ah(:, i) + facj * dx(:) + pl%ah(:, j) = pl%ah(:, j) - faci * dx(:) end do ngood = count(lgoodlevel(:)) - if (ngood > 0_I8B) then - allocate(good_idx(ngood)) - good_idx(:) = pack([(i, i = 1, nenc)], lgoodlevel(:)) - - if (isplpl) then - do concurrent (k = 1:ngood) - i = self%index1(good_idx(k)) - j = self%index2(good_idx(k)) - pl%ah(:,i) = 0.0_DP - pl%ah(:,j) = 0.0_DP - end do - else - do concurrent (k = 1_I8B:ngood) - j = self%index2(good_idx(k)) - tp%ah(:,j) = 0.0_DP - end do - end if + if (ngood == 0_I8B) return + good_idx(1:ngood) = pack([(i, i = 1, nenc)], lgoodlevel(:)) + + do k = 1, ngood + i = self%index1(good_idx(k)) + j = self%index2(good_idx(k)) + pl%vb(:,i) = pl%vb(:,i) + sgn * dt * pl%ah(:,i) + pl%vb(:,j) = pl%vb(:,j) + sgn * dt * pl%ah(:,j) + pl%ah(:,i) = 0.0_DP + pl%ah(:,j) = 0.0_DP + end do + + end if + end associate + end select + + return + end subroutine symba_kick_list_plpl + + + module subroutine symba_kick_list_pltp(self, system, dt, irec, sgn) + !! author: David A. Minton + !! + !! Kick barycentric velocities of ACTIVE test particles within SyMBA recursion. + !! + !! Adapted from David E. Kaufmann's Swifter routine: symba_kick.f90 + !! Adapted from Hal Levison's Swift routine symba5_kick.f + implicit none + ! Arguments + class(symba_list_pltp), intent(in) :: self !! SyMBA pl-tp encounter list object + class(symba_nbody_system), intent(inout) :: system !! SyMBA nbody system object + real(DP), intent(in) :: dt !! step size + integer(I4B), intent(in) :: irec !! Current recursion level + integer(I4B), intent(in) :: sgn !! sign to be applied to acceleration + ! Internals + integer(I4B) :: i, j, irm1, irecl, ngood + integer(I8B) :: k + real(DP) :: r, rr, ri, ris, rim1, r2, ir3, fac, faci + real(DP), dimension(NDIM) :: dx + logical, dimension(:), allocatable :: lgoodlevel + integer(I4B), dimension(:), allocatable :: good_idx + + if (self%nenc == 0) return + + select type(pl => system%pl) + class is (symba_pl) + select type(tp => system%tp) + class is (symba_tp) + associate(npl => pl%nbody, ntp => tp%nbody, nenc => self%nenc) + if ((npl == 0) .or. (ntp == 0)) return + pl%lmask(1:npl) = pl%status(1:npl) /= INACTIVE + tp%lmask(1:ntp) = tp%status(1:ntp) /= INACTIVE + allocate(lgoodlevel(nenc)) + + irm1 = irec - 1 + + if (sgn < 0) then + irecl = irec - 1 + else + irecl = irec + end if + + do k = 1, nenc + i = self%index1(k) + j = self%index2(k) + lgoodlevel(k) = (pl%levelg(i) >= irm1) .and. (tp%levelg(j) >= irm1) + lgoodlevel(k) = (self%status(k) == ACTIVE) .and. lgoodlevel(k) + end do + + ngood = count(lgoodlevel(:)) + + if (ngood > 0_I8B) then + allocate(good_idx(ngood)) + good_idx(:) = pack([(i, i = 1, nenc)], lgoodlevel(:)) - do k = 1, ngood - i = self%index1(good_idx(k)) - j = self%index2(good_idx(k)) - if (isplpl) then - ri = ((pl%rhill(i) + pl%rhill(j))**2) * (RHSCALE**2) * (RSHELL**(2*irecl)) - rim1 = ri * (RSHELL**2) - dx(:) = pl%rh(:,j) - pl%rh(:,i) - else - ri = ((pl%rhill(i))**2) * (RHSCALE**2) * (RSHELL**(2*irecl)) - rim1 = ri * (RSHELL**2) - dx(:) = tp%rh(:,j) - pl%rh(:,i) - end if - r2 = dot_product(dx(:), dx(:)) - if (r2 < rim1) then - fac = 0.0_DP - lgoodlevel(good_idx(k)) = .false. - cycle - end if - if (r2 < ri) then - ris = sqrt(ri) - r = sqrt(r2) - rr = (ris - r) / (ris * (1.0_DP - RSHELL)) - fac = (r2**(-1.5_DP)) * (1.0_DP - 3 * (rr**2) + 2 * (rr**3)) - else - ir3 = 1.0_DP / (r2 * sqrt(r2)) - fac = ir3 - end if - faci = fac * pl%Gmass(i) - if (isplpl) then - facj = fac * pl%Gmass(j) - pl%ah(:, i) = pl%ah(:, i) + facj * dx(:) - pl%ah(:, j) = pl%ah(:, j) - faci * dx(:) - else - tp%ah(:, j) = tp%ah(:, j) - faci * dx(:) - end if - end do - ngood = count(lgoodlevel(:)) - if (ngood == 0_I8B) return - good_idx(1:ngood) = pack([(i, i = 1, nenc)], lgoodlevel(:)) - - if (isplpl) then - do k = 1, ngood - i = self%index1(good_idx(k)) - j = self%index2(good_idx(k)) - pl%vb(:,i) = pl%vb(:,i) + sgn * dt * pl%ah(:,i) - pl%vb(:,j) = pl%vb(:,j) + sgn * dt * pl%ah(:,j) - pl%ah(:,i) = 0.0_DP - pl%ah(:,j) = 0.0_DP - end do + + do concurrent (k = 1_I8B:ngood) + j = self%index2(good_idx(k)) + tp%ah(:,j) = 0.0_DP + end do + + do k = 1, ngood + i = self%index1(good_idx(k)) + j = self%index2(good_idx(k)) + + ri = ((pl%rhill(i))**2) * (RHSCALE**2) * (RSHELL**(2*irecl)) + rim1 = ri * (RSHELL**2) + dx(:) = tp%rh(:,j) - pl%rh(:,i) + r2 = dot_product(dx(:), dx(:)) + if (r2 < rim1) then + fac = 0.0_DP + lgoodlevel(good_idx(k)) = .false. + cycle + end if + if (r2 < ri) then + ris = sqrt(ri) + r = sqrt(r2) + rr = (ris - r) / (ris * (1.0_DP - RSHELL)) + fac = (r2**(-1.5_DP)) * (1.0_DP - 3 * (rr**2) + 2 * (rr**3)) else - do k = 1, ngood - j = self%index2(good_idx(k)) - tp%vb(:,j) = tp%vb(:,j) + sgn * dt * tp%ah(:,j) - tp%ah(:,j) = 0.0_DP - end do + ir3 = 1.0_DP / (r2 * sqrt(r2)) + fac = ir3 end if - end if - end associate - end select + faci = fac * pl%Gmass(i) + + tp%ah(:, j) = tp%ah(:, j) - faci * dx(:) + end do + ngood = count(lgoodlevel(:)) + if (ngood == 0_I8B) return + good_idx(1:ngood) = pack([(i, i = 1, nenc)], lgoodlevel(:)) + + do k = 1, ngood + j = self%index2(good_idx(k)) + tp%vb(:,j) = tp%vb(:,j) + sgn * dt * tp%ah(:,j) + tp%ah(:,j) = 0.0_DP + end do + end if + end associate + end select end select return - end subroutine symba_kick_encounter + end subroutine symba_kick_list_pltp + end submodule s_symba_kick \ No newline at end of file diff --git a/src/symba/symba_setup.f90 b/src/symba/symba_setup.f90 index efd60f636..0431da5ef 100644 --- a/src/symba/symba_setup.f90 +++ b/src/symba/symba_setup.f90 @@ -70,31 +70,6 @@ module subroutine symba_setup_pl(self, n, param) end subroutine symba_setup_pl - module subroutine symba_setup_encounter_list(self, n) - !! author: David A. Minton - !! - !! A constructor that sets the number of encounters and allocates and initializes all arrays - !! - implicit none - ! Arguments - class(symba_encounter), intent(inout) :: self !! SyMBA pl-tp encounter structure - integer(I8B), intent(in) :: n !! Number of encounters to allocate space for - - call encounter_setup_list(self, n) - if (n <= 0_I8B) return - - if (allocated(self%level)) deallocate(self%level) - if (allocated(self%tcollision)) deallocate(self%tcollision) - allocate(self%level(n)) - allocate(self%tcollision(n)) - - self%level(:) = -1 - self%tcollision(:) = 0.0_DP - - return - end subroutine symba_setup_encounter_list - - module subroutine symba_setup_tp(self, n, param) !! author: David A. Minton !! diff --git a/src/symba/symba_step.f90 b/src/symba/symba_step.f90 index d51f66d89..0c4ff41f2 100644 --- a/src/symba/symba_step.f90 +++ b/src/symba/symba_step.f90 @@ -198,7 +198,7 @@ recursive module subroutine symba_step_recur_system(self, param, t, ireci) write(*, *) "SWIFTEST Warning:" write(*, *) " In symba_step_recur_system, local time step is too small" write(*, *) " Roundoff error will be important!" - call util_exit(FAILURE) + call swiftest_util_exit(FAILURE) END IF irecp = ireci + 1 if (ireci == 0) then diff --git a/src/symba/symba_util.f90 b/src/symba/symba_util.f90 index f4b143c3d..bf749a35e 100644 --- a/src/symba/symba_util.f90 +++ b/src/symba/symba_util.f90 @@ -11,32 +11,6 @@ use swiftest contains - - module subroutine symba_util_append_encounter_list(self, source, lsource_mask) - !! author: David A. Minton - !! - !! Append components from one encounter list (pl-pl or pl-tp) body object to another. - !! This method will automatically resize the destination body if it is too small - implicit none - ! Arguments - class(symba_encounter), intent(inout) :: self !! SyMBA encounter list object - class(encounter_list), intent(in) :: source !! Source object to append - logical, dimension(:), intent(in) :: lsource_mask !! Logical mask indicating which elements to append to - ! Internals - integer(I4B) :: nold, nsrc - - nold = self%nenc - nsrc = source%nenc - select type(source) - class is (symba_encounter) - call util_append(self%level, source%level, nold, nsrc, lsource_mask) - end select - call encounter_util_append_list(self, source, lsource_mask) - - return - end subroutine symba_util_append_encounter_list - - module subroutine symba_util_append_pl(self, source, lsource_mask) !! author: David A. Minton !! @@ -51,68 +25,29 @@ module subroutine symba_util_append_pl(self, source, lsource_mask) select type(source) class is (symba_pl) associate(nold => self%nbody, nsrc => source%nbody) - call util_append(self%lcollision, source%lcollision, nold, nsrc, lsource_mask) - call util_append(self%lencounter, source%lencounter, nold, nsrc, lsource_mask) - call util_append(self%lmtiny, source%lmtiny, nold, nsrc, lsource_mask) - call util_append(self%nplenc, source%nplenc, nold, nsrc, lsource_mask) - call util_append(self%ntpenc, source%ntpenc, nold, nsrc, lsource_mask) - call util_append(self%levelg, source%levelg, nold, nsrc, lsource_mask) - call util_append(self%levelm, source%levelm, nold, nsrc, lsource_mask) - call util_append(self%isperi, source%isperi, nold, nsrc, lsource_mask) - call util_append(self%peri, source%peri, nold, nsrc, lsource_mask) - call util_append(self%atp, source%atp, nold, nsrc, lsource_mask) - call util_append(self%kin, source%kin, nold, nsrc, lsource_mask) - - call util_append_pl(self, source, lsource_mask) ! Note: helio_pl does not have its own append method, so we skip back to the base class + call swiftest_util_append(self%lcollision, source%lcollision, nold, nsrc, lsource_mask) + call swiftest_util_append(self%lencounter, source%lencounter, nold, nsrc, lsource_mask) + call swiftest_util_append(self%lmtiny, source%lmtiny, nold, nsrc, lsource_mask) + call swiftest_util_append(self%nplenc, source%nplenc, nold, nsrc, lsource_mask) + call swiftest_util_append(self%ntpenc, source%ntpenc, nold, nsrc, lsource_mask) + call swiftest_util_append(self%levelg, source%levelg, nold, nsrc, lsource_mask) + call swiftest_util_append(self%levelm, source%levelm, nold, nsrc, lsource_mask) + call swiftest_util_append(self%isperi, source%isperi, nold, nsrc, lsource_mask) + call swiftest_util_append(self%peri, source%peri, nold, nsrc, lsource_mask) + call swiftest_util_append(self%atp, source%atp, nold, nsrc, lsource_mask) + call swiftest_util_append(self%kin, source%kin, nold, nsrc, lsource_mask) + + call swiftest_util_append_pl(self, source, lsource_mask) ! Note: helio_pl does not have its own append method, so we skip back to the base class end associate class default write(*,*) "Invalid object passed to the append method. Source must be of class symba_pl or its descendents!" - call util_exit(FAILURE) + call swiftest_util_exit(FAILURE) end select return end subroutine symba_util_append_pl - module subroutine symba_util_append_merger(self, source, lsource_mask) - !! author: David A. Minton - !! - !! Append components from one massive body object to another. - !! This method will automatically resize the destination body if it is too small - implicit none - ! Arguments - class(symba_merger), intent(inout) :: self !! SyMBA massive body object - class(swiftest_body), intent(in) :: source !! Source object to append - logical, dimension(:), intent(in) :: lsource_mask !! Logical mask indicating which elements to append to - ! Internals - integer(I4B), dimension(:), allocatable :: ncomp_tmp !! Temporary placeholder for ncomp incase we are appending a symba_pl object to a symba_merger - integer(I4B) :: nold, nsrc, nnew - - nold = self%nbody - nsrc = source%nbody - nnew = count(lsource_mask) - - select type(source) - class is (symba_merger) - call util_append(self%ncomp, source%ncomp, nold, nsrc, lsource_mask) - call symba_util_append_pl(self, source, lsource_mask) - class is (symba_pl) - allocate(ncomp_tmp, mold=source%id) - ncomp_tmp(:) = 0 - call util_append(self%ncomp, ncomp_tmp, nold, nsrc, lsource_mask) - call symba_util_append_pl(self, source, lsource_mask) - class default - write(*,*) "Invalid object passed to the append method. Source must be of class symba_pl or its descendents!" - call util_exit(FAILURE) - end select - - ! Save the number of appended bodies - self%ncomp(nold+1:nold+nnew) = nnew - - return - end subroutine symba_util_append_merger - - module subroutine symba_util_append_tp(self, source, lsource_mask) !! author: David A. Minton !! @@ -127,77 +62,20 @@ module subroutine symba_util_append_tp(self, source, lsource_mask) select type(source) class is (symba_tp) associate(nold => self%nbody, nsrc => source%nbody) - call util_append(self%nplenc, source%nplenc, nold, nsrc, lsource_mask) - call util_append(self%levelg, source%levelg, nold, nsrc, lsource_mask) - call util_append(self%levelm, source%levelm, nold, nsrc, lsource_mask) + call swiftest_util_append(self%levelg, source%levelg, nold, nsrc, lsource_mask) + call swiftest_util_append(self%levelm, source%levelm, nold, nsrc, lsource_mask) - call util_append_tp(self, source, lsource_mask) ! Note: helio_tp does not have its own append method, so we skip back to the base class + call swiftest_util_append_tp(self, source, lsource_mask) ! Note: helio_tp does not have its own append method, so we skip back to the base class end associate class default write(*,*) "Invalid object passed to the append method. Source must be of class symba_tp or its descendents!" - call util_exit(FAILURE) + call swiftest_util_exit(FAILURE) end select return end subroutine symba_util_append_tp - module subroutine symba_util_copy_encounter_list(self, source) - !! author: David A. Minton - !! - !! Copies elements from the source encounter list into self. - implicit none - ! Arguments - class(symba_encounter), intent(inout) :: self !! Encounter list - class(encounter_list), intent(in) :: source !! Source object to copy into - - select type(source) - class is (symba_encounter) - associate(n => source%nenc) - self%level(1:n) = source%level(1:n) - self%tcollision(1:n) = source%tcollision(1:n) - end associate - end select - - call encounter_util_copy_list(self, source) - - return - end subroutine symba_util_copy_encounter_list - - - module subroutine symba_util_dealloc_encounter_list(self) - !! author: David A. Minton - !! - !! Deallocates all allocatabale arrays - implicit none - ! Argumentse - class(symba_encounter), intent(inout) :: self !! SyMBA encounter list - - if (allocated(self%level)) deallocate(self%level) - if (allocated(self%tcollision)) deallocate(self%tcollision) - - call self%encounter_list%dealloc() - - return - end subroutine symba_util_dealloc_encounter_list - - - module subroutine symba_util_dealloc_merger(self) - !! author: David A. Minton - !! - !! Deallocates all allocatabale arrays - implicit none - ! Arguments - class(symba_merger), intent(inout) :: self !! SyMBA body merger object - - if (allocated(self%ncomp)) deallocate(self%ncomp) - - call self%symba_pl%dealloc() - - return - end subroutine symba_util_dealloc_merger - - module subroutine symba_util_dealloc_pl(self) !! author: David A. Minton !! @@ -252,22 +130,22 @@ module subroutine symba_util_fill_pl(self, inserts, lfill_list) associate(keeps => self) select type(inserts) class is (symba_pl) - call util_fill(keeps%lcollision, inserts%lcollision, lfill_list) - call util_fill(keeps%lencounter, inserts%lencounter, lfill_list) - call util_fill(keeps%lmtiny, inserts%lmtiny, lfill_list) - call util_fill(keeps%nplenc, inserts%nplenc, lfill_list) - call util_fill(keeps%ntpenc, inserts%ntpenc, lfill_list) - call util_fill(keeps%levelg, inserts%levelg, lfill_list) - call util_fill(keeps%levelm, inserts%levelm, lfill_list) - call util_fill(keeps%isperi, inserts%isperi, lfill_list) - call util_fill(keeps%peri, inserts%peri, lfill_list) - call util_fill(keeps%atp, inserts%atp, lfill_list) - call util_fill(keeps%kin, inserts%kin, lfill_list) + call swiftest_util_fill(keeps%lcollision, inserts%lcollision, lfill_list) + call swiftest_util_fill(keeps%lencounter, inserts%lencounter, lfill_list) + call swiftest_util_fill(keeps%lmtiny, inserts%lmtiny, lfill_list) + call swiftest_util_fill(keeps%nplenc, inserts%nplenc, lfill_list) + call swiftest_util_fill(keeps%ntpenc, inserts%ntpenc, lfill_list) + call swiftest_util_fill(keeps%levelg, inserts%levelg, lfill_list) + call swiftest_util_fill(keeps%levelm, inserts%levelm, lfill_list) + call swiftest_util_fill(keeps%isperi, inserts%isperi, lfill_list) + call swiftest_util_fill(keeps%peri, inserts%peri, lfill_list) + call swiftest_util_fill(keeps%atp, inserts%atp, lfill_list) + call swiftest_util_fill(keeps%kin, inserts%kin, lfill_list) - call util_fill_pl(keeps, inserts, lfill_list) ! Note: helio_pl does not have its own fill method, so we skip back to the base class + call swiftest_util_fill_pl(keeps, inserts, lfill_list) ! Note: helio_pl does not have its own fill method, so we skip back to the base class class default write(*,*) "Invalid object passed to the fill method. Source must be of class symba_pl or its descendents!" - call util_exit(FAILURE) + call swiftest_util_exit(FAILURE) end select end associate @@ -290,14 +168,14 @@ module subroutine symba_util_fill_tp(self, inserts, lfill_list) associate(keeps => self) select type(inserts) class is (symba_tp) - call util_fill(keeps%nplenc, inserts%nplenc, lfill_list) - call util_fill(keeps%levelg, inserts%levelg, lfill_list) - call util_fill(keeps%levelm, inserts%levelm, lfill_list) + call swiftest_util_fill(keeps%nplenc, inserts%nplenc, lfill_list) + call swiftest_util_fill(keeps%levelg, inserts%levelg, lfill_list) + call swiftest_util_fill(keeps%levelm, inserts%levelm, lfill_list) - call util_fill_tp(keeps, inserts, lfill_list) ! Note: helio_tp does not have its own fill method, so we skip back to the base class + call swiftest_util_fill_tp(keeps, inserts, lfill_list) ! Note: helio_tp does not have its own fill method, so we skip back to the base class class default write(*,*) "Invalid object passed to the fill method. Source must be of class symba_tp or its descendents!" - call util_exit(FAILURE) + call swiftest_util_exit(FAILURE) end select end associate @@ -305,67 +183,30 @@ module subroutine symba_util_fill_tp(self, inserts, lfill_list) end subroutine symba_util_fill_tp - module subroutine symba_util_flatten_eucl_plpl(self, param) - !! author: Jacob R. Elliott and David A. Minton - !! - !! Turns i,j indices into k index for use in the Euclidean distance matrix. This also sets the lmtiny flag and computes the - !! number of interactions that excludes semi-interacting bodies with each other (Gmass < GMTINY). - !! This method will also sort the bodies in descending order by Mass - !! - !! Reference: - !! - !! Mélodie Angeletti, Jean-Marie Bonny, Jonas Koko. Parallel Euclidean distance matrix computation on big datasets *. - !! 2019. hal-0204751 - implicit none - ! Arguments - class(symba_pl), intent(inout) :: self !! SyMBA massive body object - class(swiftest_parameters), intent(inout) :: param !! Current run configuration parameters - ! Internals - integer(I8B) :: npl, nplm - - associate(pl => self, nplplm => self%nplplm) - npl = int(self%nbody, kind=I8B) - select type(param) - class is (swiftest_parameters) - pl%lmtiny(1:npl) = pl%Gmass(1:npl) < param%GMTINY - end select - nplm = count(.not. pl%lmtiny(1:npl)) - pl%nplm = int(nplm, kind=I4B) - nplplm = nplm * npl - nplm * (nplm + 1_I8B) / 2_I8B ! number of entries in a strict lower triangle, npl x npl, minus first column including only mutually interacting bodies - - call util_flatten_eucl_plpl(pl, param) - end associate - - return - end subroutine symba_util_flatten_eucl_plpl - - - module subroutine symba_util_final_encounter_list(self) + module subroutine symba_util_final_list_plpl(self) !! author: David A. Minton !! - !! Finalize the SyMBA encounter list object - deallocates all allocatables + !! Finalize the pl-tp list - deallocates all allocatables implicit none - ! Argument - type(symba_encounter), intent(inout) :: self !! SyMBA encounter list object + type(symba_list_plpl), intent(inout) :: self !! SyMBA encounter list object call self%dealloc() return - end subroutine symba_util_final_encounter_list + end subroutine symba_util_final_list_plpl - module subroutine symba_util_final_merger(self) + module subroutine symba_util_final_list_pltp(self) !! author: David A. Minton !! - !! Finalize the SyMBA merger object - deallocates all allocatables + !! Finalize the pl-tp list - deallocates all allocatables implicit none - ! Argument - type(symba_merger), intent(inout) :: self !! SyMBA merger object + type(symba_list_pltp), intent(inout) :: self !! SyMBA encounter list object call self%dealloc() return - end subroutine symba_util_final_merger + end subroutine symba_util_final_list_pltp module subroutine symba_util_final_pl(self) @@ -415,6 +256,41 @@ module subroutine symba_util_final_tp(self) end subroutine symba_util_final_tp + module subroutine symba_util_flatten_eucl_plpl(self, param) + !! author: Jacob R. Elliott and David A. Minton + !! + !! Turns i,j indices into k index for use in the Euclidean distance matrix. This also sets the lmtiny flag and computes the + !! number of interactions that excludes semi-interacting bodies with each other (Gmass < GMTINY). + !! This method will also sort the bodies in descending order by Mass + !! + !! Reference: + !! + !! Mélodie Angeletti, Jean-Marie Bonny, Jonas Koko. Parallel Euclidean distance matrix computation on big datasets *. + !! 2019. hal-0204751 + implicit none + ! Arguments + class(symba_pl), intent(inout) :: self !! SyMBA massive body object + class(swiftest_parameters), intent(inout) :: param !! Current run configuration parameters + ! Internals + integer(I8B) :: npl, nplm + + associate(pl => self, nplplm => self%nplplm) + npl = int(self%nbody, kind=I8B) + select type(param) + class is (swiftest_parameters) + pl%lmtiny(1:npl) = pl%Gmass(1:npl) < param%GMTINY + end select + nplm = count(.not. pl%lmtiny(1:npl)) + pl%nplm = int(nplm, kind=I4B) + nplplm = nplm * npl - nplm * (nplm + 1_I8B) / 2_I8B ! number of entries in a strict lower triangle, npl x npl, minus first column including only mutually interacting bodies + + call swiftest_util_flatten_eucl_plpl(pl, param) + end associate + + return + end subroutine symba_util_flatten_eucl_plpl + + module subroutine symba_util_resize_pl(self, nnew) !! author: David A. Minton !! @@ -424,13 +300,13 @@ module subroutine symba_util_resize_pl(self, nnew) class(symba_pl), intent(inout) :: self !! SyMBA massive body object integer(I4B), intent(in) :: nnew !! New size neded - call util_resize(self%levelg, nnew) - call util_resize(self%levelm, nnew) - call util_resize(self%isperi, nnew) - call util_resize(self%peri, nnew) - call util_resize(self%atp, nnew) + call swiftest_util_resize(self%levelg, nnew) + call swiftest_util_resize(self%levelm, nnew) + call swiftest_util_resize(self%isperi, nnew) + call swiftest_util_resize(self%peri, nnew) + call swiftest_util_resize(self%atp, nnew) - call util_resize_pl(self, nnew) + call swiftest_util_resize_pl(self, nnew) return end subroutine symba_util_resize_pl @@ -444,10 +320,10 @@ module subroutine symba_util_resize_tp(self, nnew) class(symba_tp), intent(inout) :: self !! SyMBA test particle object integer(I4B), intent(in) :: nnew !! New size neded - call util_resize(self%levelg, nnew) - call util_resize(self%levelm, nnew) + call swiftest_util_resize(self%levelg, nnew) + call swiftest_util_resize(self%levelm, nnew) - call util_resize_tp(self, nnew) + call swiftest_util_resize_tp(self, nnew) return end subroutine symba_util_resize_tp @@ -503,21 +379,21 @@ module subroutine symba_util_sort_pl(self, sortby, ascending) associate(pl => self, npl => self%nbody) select case(sortby) case("nplenc") - call util_sort(direction * pl%nplenc(1:npl), ind) + call swiftest_util_sort(direction * pl%nplenc(1:npl), ind) case("ntpenc") - call util_sort(direction * pl%ntpenc(1:npl), ind) + call swiftest_util_sort(direction * pl%ntpenc(1:npl), ind) case("levelg") - call util_sort(direction * pl%levelg(1:npl), ind) + call swiftest_util_sort(direction * pl%levelg(1:npl), ind) case("levelm") - call util_sort(direction * pl%levelm(1:npl), ind) + call swiftest_util_sort(direction * pl%levelm(1:npl), ind) case("peri") - call util_sort(direction * pl%peri(1:npl), ind) + call swiftest_util_sort(direction * pl%peri(1:npl), ind) case("atp") - call util_sort(direction * pl%atp(1:npl), ind) + call swiftest_util_sort(direction * pl%atp(1:npl), ind) case("lcollision", "lencounter", "lmtiny", "nplm", "nplplm", "kin", "info") write(*,*) 'Cannot sort by ' // trim(adjustl(sortby)) // '. Component not sortable!' case default ! Look for components in the parent class - call util_sort_pl(pl, sortby, ascending) + call swiftest_util_sort_pl(pl, sortby, ascending) return end select @@ -553,13 +429,13 @@ module subroutine symba_util_sort_tp(self, sortby, ascending) associate(tp => self, ntp => self%nbody) select case(sortby) case("nplenc") - call util_sort(direction * tp%nplenc(1:ntp), ind) + call swiftest_util_sort(direction * tp%nplenc(1:ntp), ind) case("levelg") - call util_sort(direction * tp%levelg(1:ntp), ind) + call swiftest_util_sort(direction * tp%levelg(1:ntp), ind) case("levelm") - call util_sort(direction * tp%levelm(1:ntp), ind) + call swiftest_util_sort(direction * tp%levelm(1:ntp), ind) case default ! Look for components in the parent class - call util_sort_tp(tp, sortby, ascending) + call swiftest_util_sort_tp(tp, sortby, ascending) return end select @@ -581,9 +457,9 @@ module subroutine symba_util_sort_rearrange_pl(self, ind) integer(I4B), dimension(:), intent(in) :: ind !! Index array used to restructure the body (should contain all 1:n index values in the desired order) associate(pl => self, npl => self%nbody) - call util_sort_rearrange(pl%levelg, ind, npl) - call util_sort_rearrange(pl%levelm, ind, npl) - call util_sort_rearrange_pl(pl,ind) + call swiftest_util_sort_rearrange(pl%levelg, ind, npl) + call swiftest_util_sort_rearrange(pl%levelm, ind, npl) + call swiftest_util_sort_rearrange_pl(pl,ind) end associate return @@ -601,19 +477,17 @@ module subroutine symba_util_sort_rearrange_tp(self, ind) integer(I4B), dimension(:), intent(in) :: ind !! Index array used to restructure the body (should contain all 1:n index values in the desired order) associate(tp => self, ntp => self%nbody) - call util_sort_rearrange(tp%nplenc, ind, ntp) - call util_sort_rearrange(tp%levelg, ind, ntp) - call util_sort_rearrange(tp%levelm, ind, ntp) + call swiftest_util_sort_rearrange(tp%nplenc, ind, ntp) + call swiftest_util_sort_rearrange(tp%levelg, ind, ntp) + call swiftest_util_sort_rearrange(tp%levelm, ind, ntp) - call util_sort_rearrange_tp(tp,ind) + call swiftest_util_sort_rearrange_tp(tp,ind) end associate return end subroutine symba_util_sort_rearrange_tp - - module subroutine symba_util_spill_pl(self, discards, lspill_list, ldestructive) !! author: David A. Minton !! @@ -631,22 +505,22 @@ module subroutine symba_util_spill_pl(self, discards, lspill_list, ldestructive) associate(keeps => self) select type(discards) class is (symba_pl) - call util_spill(keeps%lcollision, discards%lcollision, lspill_list, ldestructive) - call util_spill(keeps%lencounter, discards%lencounter, lspill_list, ldestructive) - call util_spill(keeps%lmtiny, discards%lmtiny, lspill_list, ldestructive) - call util_spill(keeps%nplenc, discards%nplenc, lspill_list, ldestructive) - call util_spill(keeps%ntpenc, discards%ntpenc, lspill_list, ldestructive) - call util_spill(keeps%levelg, discards%levelg, lspill_list, ldestructive) - call util_spill(keeps%levelm, discards%levelm, lspill_list, ldestructive) - call util_spill(keeps%isperi, discards%isperi, lspill_list, ldestructive) - call util_spill(keeps%peri, discards%peri, lspill_list, ldestructive) - call util_spill(keeps%atp, discards%atp, lspill_list, ldestructive) - call util_spill(keeps%kin, discards%kin, lspill_list, ldestructive) - - call util_spill_pl(keeps, discards, lspill_list, ldestructive) + call swiftest_util_spill(keeps%lcollision, discards%lcollision, lspill_list, ldestructive) + call swiftest_util_spill(keeps%lencounter, discards%lencounter, lspill_list, ldestructive) + call swiftest_util_spill(keeps%lmtiny, discards%lmtiny, lspill_list, ldestructive) + call swiftest_util_spill(keeps%nplenc, discards%nplenc, lspill_list, ldestructive) + call swiftest_util_spill(keeps%ntpenc, discards%ntpenc, lspill_list, ldestructive) + call swiftest_util_spill(keeps%levelg, discards%levelg, lspill_list, ldestructive) + call swiftest_util_spill(keeps%levelm, discards%levelm, lspill_list, ldestructive) + call swiftest_util_spill(keeps%isperi, discards%isperi, lspill_list, ldestructive) + call swiftest_util_spill(keeps%peri, discards%peri, lspill_list, ldestructive) + call swiftest_util_spill(keeps%atp, discards%atp, lspill_list, ldestructive) + call swiftest_util_spill(keeps%kin, discards%kin, lspill_list, ldestructive) + + call swiftest_util_spill_pl(keeps, discards, lspill_list, ldestructive) class default write(*,*) "Invalid object passed to the spill method. Source must be of class symba_pl or its descendents!" - call util_exit(FAILURE) + call swiftest_util_exit(FAILURE) end select end associate @@ -654,32 +528,6 @@ module subroutine symba_util_spill_pl(self, discards, lspill_list, ldestructive) end subroutine symba_util_spill_pl - module subroutine symba_util_spill_encounter_list(self, discards, lspill_list, ldestructive) - !! author: David A. Minton - !! - !! Move spilled (discarded) SyMBA encounter structure from active list to discard list - !! Note: Because the plpl_encounter currently does not contain any additional variable components, this method can recieve it as an input as well. - implicit none - ! Arguments - class(symba_encounter), intent(inout) :: self !! SyMBA pl-tp encounter list - class(encounter_list), intent(inout) :: discards !! Discarded object - logical, dimension(:), intent(in) :: lspill_list !! Logical array of bodies to spill into the discards - logical, intent(in) :: ldestructive !! Logical flag indicating whether or not this operation should alter body by removing the discard list - - associate(keeps => self) - select type(discards) - class is (symba_encounter) - call util_spill(keeps%level, discards%level, lspill_list, ldestructive) - call encounter_util_spill_list(keeps, discards, lspill_list, ldestructive) - class default - write(*,*) "Invalid object passed to the spill method. Source must be of class symba_encounter or its descendents!" - call util_exit(FAILURE) - end select - end associate - - return - end subroutine symba_util_spill_encounter_list - module subroutine symba_util_spill_tp(self, discards, lspill_list, ldestructive) !! author: David A. Minton @@ -698,14 +546,14 @@ module subroutine symba_util_spill_tp(self, discards, lspill_list, ldestructive) associate(keeps => self) select type(discards) class is (symba_tp) - call util_spill(keeps%nplenc, discards%nplenc, lspill_list, ldestructive) - call util_spill(keeps%levelg, discards%levelg, lspill_list, ldestructive) - call util_spill(keeps%levelm, discards%levelm, lspill_list, ldestructive) + call swiftest_util_spill(keeps%nplenc, discards%nplenc, lspill_list, ldestructive) + call swiftest_util_spill(keeps%levelg, discards%levelg, lspill_list, ldestructive) + call swiftest_util_spill(keeps%levelm, discards%levelm, lspill_list, ldestructive) - call util_spill_tp(keeps, discards, lspill_list, ldestructive) + call swiftest_util_spill_tp(keeps, discards, lspill_list, ldestructive) class default write(*,*) "Invalid object passed to the spill method. Source must be of class symba_tp or its descendents!" - call util_exit(FAILURE) + call swiftest_util_exit(FAILURE) end select end associate diff --git a/src/tides/tides_spin_step.f90 b/src/tides/tides_spin_step.f90 index 3fadfc704..fd74e4a21 100644 --- a/src/tides/tides_spin_step.f90 +++ b/src/tides/tides_spin_step.f90 @@ -51,7 +51,7 @@ module subroutine tides_step_spin_system(self, param, t, dt) rot0 = [pack(pl%rot(:,1:npl),.true.), pack(cb%rot(:),.true.)] ! Use this space call the ode_solver, passing tides_spin_derivs as the function: subdt = dt / 20._DP - !rot1(:) = util_solve_rkf45(lambda_obj(tides_spin_derivs, subdt, pl%rbeg, pl%rend), rot0, dt, subdt tol) + !rot1(:) = swiftest_util_solve_rkf45(lambda_obj(tides_spin_derivs, subdt, pl%rbeg, pl%rend), rot0, dt, subdt tol) ! Recover with unpack !pl%rot(:,1:npl) = unpack(rot1... !cb%rot(:) = unpack(rot1... diff --git a/src/whm/whm_drift.f90 b/src/whm/whm_drift.f90 index 31d041505..5565ea8e8 100644 --- a/src/whm/whm_drift.f90 +++ b/src/whm/whm_drift.f90 @@ -48,7 +48,7 @@ module subroutine whm_drift_pl(self, system, param, dt) WRITE(*, *) " STOPPING " end if end do - call util_exit(FAILURE) + call swiftest_util_exit(FAILURE) end if end associate diff --git a/src/whm/whm_setup.f90 b/src/whm/whm_setup.f90 index ae54fa9e8..c6be76904 100644 --- a/src/whm/whm_setup.f90 +++ b/src/whm/whm_setup.f90 @@ -56,7 +56,7 @@ module subroutine whm_util_set_mu_eta_pl(self, cb) associate(pl => self, npl => self%nbody) if (npl == 0) return - call util_set_mu_pl(pl, cb) + call swiftest_util_set_mu_pl(pl, cb) pl%eta(1) = cb%Gmass + pl%Gmass(1) pl%muj(1) = pl%eta(1) do i = 2, npl @@ -81,7 +81,7 @@ module subroutine whm_setup_initialize_system(self, param) call setup_initialize_system(self, param) ! First we need to make sure that the massive bodies are sorted by heliocentric distance before computing jacobies - call util_set_ir3h(self%pl) + call swiftest_util_set_ir3h(self%pl) call self%pl%sort("ir3h", ascending=.false.) call self%pl%flatten(param) diff --git a/src/whm/whm_util.f90 b/src/whm/whm_util.f90 index f6a16a44c..6d49c818f 100644 --- a/src/whm/whm_util.f90 +++ b/src/whm/whm_util.f90 @@ -25,17 +25,17 @@ module subroutine whm_util_append_pl(self, source, lsource_mask) select type(source) class is (whm_pl) associate(nold => self%nbody, nsrc => source%nbody) - call util_append(self%eta, source%eta, nold, nsrc, lsource_mask) - call util_append(self%muj, source%muj, nold, nsrc, lsource_mask) - call util_append(self%ir3j, source%ir3j, nold, nsrc, lsource_mask) - call util_append(self%xj, source%xj, nold, nsrc, lsource_mask) - call util_append(self%vj, source%vj, nold, nsrc, lsource_mask) + call swiftest_util_append(self%eta, source%eta, nold, nsrc, lsource_mask) + call swiftest_util_append(self%muj, source%muj, nold, nsrc, lsource_mask) + call swiftest_util_append(self%ir3j, source%ir3j, nold, nsrc, lsource_mask) + call swiftest_util_append(self%xj, source%xj, nold, nsrc, lsource_mask) + call swiftest_util_append(self%vj, source%vj, nold, nsrc, lsource_mask) - call util_append_pl(self, source, lsource_mask) + call swiftest_util_append_pl(self, source, lsource_mask) end associate class default write(*,*) "Invalid object passed to the append method. Source must be of class whm_pl or its descendents" - call util_exit(FAILURE) + call swiftest_util_exit(FAILURE) end select return @@ -56,7 +56,7 @@ module subroutine whm_util_dealloc_pl(self) if (allocated(self%vj)) deallocate(self%vj) if (allocated(self%ir3j)) deallocate(self%ir3j) - call util_dealloc_pl(self) + call swiftest_util_dealloc_pl(self) return end subroutine whm_util_dealloc_pl @@ -78,16 +78,16 @@ module subroutine whm_util_fill_pl(self, inserts, lfill_list) associate(keeps => self) select type(inserts) class is (whm_pl) - call util_fill(keeps%eta, inserts%eta, lfill_list) - call util_fill(keeps%muj, inserts%muj, lfill_list) - call util_fill(keeps%ir3j, inserts%ir3j, lfill_list) - call util_fill(keeps%xj, inserts%xj, lfill_list) - call util_fill(keeps%vj, inserts%vj, lfill_list) + call swiftest_util_fill(keeps%eta, inserts%eta, lfill_list) + call swiftest_util_fill(keeps%muj, inserts%muj, lfill_list) + call swiftest_util_fill(keeps%ir3j, inserts%ir3j, lfill_list) + call swiftest_util_fill(keeps%xj, inserts%xj, lfill_list) + call swiftest_util_fill(keeps%vj, inserts%vj, lfill_list) - call util_fill_pl(keeps, inserts, lfill_list) + call swiftest_util_fill_pl(keeps, inserts, lfill_list) class default write(*,*) "Invalid object passed to the fill method. Inserts must be of class whm_pl or its descendents!" - call util_exit(FAILURE) + call swiftest_util_exit(FAILURE) end select end associate @@ -117,7 +117,7 @@ module subroutine whm_util_final_system(self) ! Arguments type(whm_nbody_system), intent(inout) :: self !! WHM nbody system object - call util_final_system(self) + call swiftest_util_final_system(self) return end subroutine whm_util_final_system @@ -146,13 +146,13 @@ module subroutine whm_util_resize_pl(self, nnew) class(whm_pl), intent(inout) :: self !! WHM massive body object integer(I4B), intent(in) :: nnew !! New size neded - call util_resize(self%eta, nnew) - call util_resize(self%xj, nnew) - call util_resize(self%vj, nnew) - call util_resize(self%muj, nnew) - call util_resize(self%ir3j, nnew) + call swiftest_util_resize(self%eta, nnew) + call swiftest_util_resize(self%xj, nnew) + call swiftest_util_resize(self%vj, nnew) + call swiftest_util_resize(self%muj, nnew) + call swiftest_util_resize(self%ir3j, nnew) - call util_resize_pl(self, nnew) + call swiftest_util_resize_pl(self, nnew) return end subroutine whm_util_resize_pl @@ -209,15 +209,15 @@ module subroutine whm_util_sort_pl(self, sortby, ascending) associate(pl => self, npl => self%nbody) select case(sortby) case("eta") - call util_sort(direction * pl%eta(1:npl), ind) + call swiftest_util_sort(direction * pl%eta(1:npl), ind) case("muj") - call util_sort(direction * pl%muj(1:npl), ind) + call swiftest_util_sort(direction * pl%muj(1:npl), ind) case("ir3j") - call util_sort(direction * pl%ir3j(1:npl), ind) + call swiftest_util_sort(direction * pl%ir3j(1:npl), ind) case("xj", "vj") write(*,*) 'Cannot sort by ' // trim(adjustl(sortby)) // '. Component not sortable!' case default - call util_sort_pl(pl, sortby, ascending) + call swiftest_util_sort_pl(pl, sortby, ascending) return end select @@ -241,13 +241,13 @@ module subroutine whm_util_sort_rearrange_pl(self, ind) if (self%nbody == 0) return associate(pl => self, npl => self%nbody) - call util_sort_rearrange(pl%eta, ind, npl) - call util_sort_rearrange(pl%xj, ind, npl) - call util_sort_rearrange(pl%vj, ind, npl) - call util_sort_rearrange(pl%muj, ind, npl) - call util_sort_rearrange(pl%ir3j, ind, npl) + call swiftest_util_sort_rearrange(pl%eta, ind, npl) + call swiftest_util_sort_rearrange(pl%xj, ind, npl) + call swiftest_util_sort_rearrange(pl%vj, ind, npl) + call swiftest_util_sort_rearrange(pl%muj, ind, npl) + call swiftest_util_sort_rearrange(pl%ir3j, ind, npl) - call util_sort_rearrange_pl(pl,ind) + call swiftest_util_sort_rearrange_pl(pl,ind) end associate return @@ -270,16 +270,16 @@ module subroutine whm_util_spill_pl(self, discards, lspill_list, ldestructive) associate(keeps => self) select type(discards) class is (whm_pl) - call util_spill(keeps%eta, discards%eta, lspill_list, ldestructive) - call util_spill(keeps%muj, discards%muj, lspill_list, ldestructive) - call util_spill(keeps%ir3j, discards%ir3j, lspill_list, ldestructive) - call util_spill(keeps%xj, discards%xj, lspill_list, ldestructive) - call util_spill(keeps%vj, discards%vj, lspill_list, ldestructive) + call swiftest_util_spill(keeps%eta, discards%eta, lspill_list, ldestructive) + call swiftest_util_spill(keeps%muj, discards%muj, lspill_list, ldestructive) + call swiftest_util_spill(keeps%ir3j, discards%ir3j, lspill_list, ldestructive) + call swiftest_util_spill(keeps%xj, discards%xj, lspill_list, ldestructive) + call swiftest_util_spill(keeps%vj, discards%vj, lspill_list, ldestructive) - call util_spill_pl(keeps, discards, lspill_list, ldestructive) + call swiftest_util_spill_pl(keeps, discards, lspill_list, ldestructive) class default write(*,*) "Invalid object passed to the spill method. Source must be of class whm_pl or its descendents!" - call util_exit(FAILURE) + call swiftest_util_exit(FAILURE) end select end associate