diff --git a/Makefile.Defines b/Makefile.Defines index 291f2c604..1346f4b09 100644 --- a/Makefile.Defines +++ b/Makefile.Defines @@ -65,13 +65,13 @@ GPAR = -fopenmp -ftree-parallelize-loops=4 GMEM = -fsanitize=undefined -fsanitize=address -fsanitize=leak GWARNINGS = -Wall -Warray-bounds -Wimplicit-interface -Wextra -Warray-temporaries -FFLAGS = $(IDEBUG) $(HEAPARR) +#FFLAGS = $(IDEBUG) $(HEAPARR) #FFLAGS = -init=snan,arrays -no-wrap-margin -O3 $(STRICTREAL) $(SIMDVEC) $(PAR) -FORTRAN = ifort +#FORTRAN = ifort #AR = xiar -#FORTRAN = gfortran -#FFLAGS = -ffree-line-length-none $(GDEBUG) $(GMEM) +FORTRAN = gfortran +FFLAGS = -ffree-line-length-none $(GDEBUG) #$(GMEM) AR = ar # DO NOT include in CFLAGS the "-c" option to compile object only diff --git a/src/discard/discard.f90 b/src/discard/discard.f90 index 95d4ca4b4..292e52c38 100644 --- a/src/discard/discard.f90 +++ b/src/discard/discard.f90 @@ -8,13 +8,19 @@ module subroutine discard_system(self, param) !! Calls the discard methods for each body class and then the write method if any discards were detected !! implicit none + ! Arguments class(swiftest_nbody_system), intent(inout) :: self !! Swiftest system object class(swiftest_parameters), intent(in) :: param !! Current run configuration parameters + ! Internals + logical :: lany_discards associate(system => self, tp => self%tp, pl => self%pl) - call tp%discard(system, param) + lany_discards = .false. call pl%discard(system, param) - if (any(tp%ldiscard(:) .or. any(pl%ldiscard(:)))) call system%write_discard(param) + call tp%discard(system, param) + if (tp%nbody > 0) lany_discards = lany_discards .or. any(tp%ldiscard(:)) + if (pl%nbody > 0) lany_discards = lany_discards .or. any(pl%ldiscard(:)) + if (lany_discards) call system%write_discard(param) end associate return @@ -31,6 +37,7 @@ module subroutine discard_pl(self, system, param) class(swiftest_pl), intent(inout) :: self !! Swiftest massive body object class(swiftest_nbody_system), intent(inout) :: system !! Swiftest nbody system object class(swiftest_parameters), intent(in) :: param !! Current run configuration parameter + self%ldiscard(:) = .false. return diff --git a/src/modules/symba_classes.f90 b/src/modules/symba_classes.f90 index ff2f08dc5..fb5f6ec06 100644 --- a/src/modules/symba_classes.f90 +++ b/src/modules/symba_classes.f90 @@ -387,9 +387,10 @@ module subroutine symba_step_interp_system(self, param, t, dt) real(DP), intent(in) :: dt !! Current stepsize end subroutine symba_step_interp_system - module subroutine symba_step_set_recur_levels_system(self) + module subroutine symba_step_set_recur_levels_system(self, ireci) implicit none - class(symba_nbody_system), intent(inout) :: self !! SyMBA nbody system objec + class(symba_nbody_system), intent(inout) :: self !! SyMBA nbody system objec + integer(I4B), intent(in) :: ireci !! Input recursion level end subroutine symba_step_set_recur_levels_system module recursive subroutine symba_step_recur_system(self, param, t, ireci) diff --git a/src/symba/symba_encounter_check.f90 b/src/symba/symba_encounter_check.f90 index ce3855e02..1a9d8c68f 100644 --- a/src/symba/symba_encounter_check.f90 +++ b/src/symba/symba_encounter_check.f90 @@ -17,6 +17,7 @@ module function symba_encounter_check_pl(self, system, dt, irec) result(lany_enc logical :: lany_encounter !! Returns true if there is at least one close encounter ! Internals integer(I8B) :: k + integer(I4B) :: nenc real(DP), dimension(NDIM) :: xr, vr logical, dimension(:), allocatable :: lencounter, loc_lvdotr @@ -34,17 +35,20 @@ module function symba_encounter_check_pl(self, system, dt, irec) result(lany_enc end associate end do - lany_encounter = any(lencounter(:)) + nenc = count(lencounter(:)) + lany_encounter = nenc > 0 if (lany_encounter) then - associate(plplenc_list => system%plplenc_list, nenc => system%plplenc_list%nenc) - call plplenc_list%resize(count(lencounter(:))) - plplenc_list%status(1:nenc) = ACTIVE - plplenc_list%level(1:nenc) = irec + associate(plplenc_list => system%plplenc_list) + call plplenc_list%resize(nenc) plplenc_list%lvdotr(1:nenc) = pack(loc_lvdotr(:), lencounter(:)) plplenc_list%index1(1:nenc) = pack(pl%k_plpl(1,:), lencounter(:)) plplenc_list%index2(1:nenc) = pack(pl%k_plpl(2,:), lencounter(:)) - pl%lencounter(plplenc_list%index1(1:nenc)) = .true. - pl%lencounter(plplenc_list%index2(1:nenc)) = .true. + do k = 1, nenc + plplenc_list%status(k) = ACTIVE + plplenc_list%level(k) = irec + pl%lencounter(plplenc_list%index1(k)) = .true. + pl%lencounter(plplenc_list%index2(k)) = .true. + end do end associate end if end associate @@ -142,7 +146,7 @@ module function symba_encounter_check_tp(self, system, dt, irec) result(lany_enc logical :: lany_encounter !! Returns true if there is at least one close encounter ! Internals real(DP) :: r2crit, vdotr, r2, v2, tmin, r2min, term2 - integer(I4B) :: i, j + integer(I4B) :: i, j, k,nenc real(DP), dimension(NDIM) :: xr, vr logical, dimension(:,:), allocatable :: lencounter, loc_lvdotr @@ -160,10 +164,11 @@ module function symba_encounter_check_tp(self, system, dt, irec) result(lany_enc end do end do - lany_encounter = any(lencounter(:,:)) + nenc = count(lencounter(:,:)) + lany_encounter = nenc > 0 if (lany_encounter) then - associate(pltpenc_list => system%pltpenc_list, nenc => system%pltpenc_list%nenc) - call pltpenc_list%resize(count(lencounter(:,:))) + associate(pltpenc_list => system%pltpenc_list) + call pltpenc_list%resize(nenc) pltpenc_list%status(1:nenc) = ACTIVE pltpenc_list%level(1:nenc) = irec pltpenc_list%lvdotr(1:nenc) = pack(loc_lvdotr(:,:), lencounter(:,:)) @@ -172,7 +177,9 @@ module function symba_encounter_check_tp(self, system, dt, irec) result(lany_enc select type(pl) class is (symba_pl) pl%lencounter(:) = .false. - pl%lencounter(pltpenc_list%index1(1:nenc)) = .true. + do k = 1, nenc + pl%lencounter(pltpenc_list%index1(k)) = .true. + end do end select end associate end if diff --git a/src/symba/symba_step.f90 b/src/symba/symba_step.f90 index 73319033c..c34f9fb0c 100644 --- a/src/symba/symba_step.f90 +++ b/src/symba/symba_step.f90 @@ -94,14 +94,15 @@ module subroutine symba_step_interp_system(self, param, t, dt) end subroutine symba_step_interp_system - module subroutine symba_step_set_recur_levels_system(self) + module subroutine symba_step_set_recur_levels_system(self, ireci) !! author: David A. Minton !! !! Resets pl, tp,and encounter structures at the start of a new step !! implicit none ! Arguments - class(symba_nbody_system), intent(inout) :: self !! SyMBA nbody system object + class(symba_nbody_system), intent(inout) :: self !! SyMBA nbody system object + integer(I4B), intent(in) :: ireci !! Input recursion level ! Internals integer(I4B) :: i, irecp @@ -115,20 +116,21 @@ module subroutine symba_step_set_recur_levels_system(self) plind3 => pltpenc_list%index1(1:pltpenc_list%nenc), & tpind => pltpenc_list%index2(1:pltpenc_list%nenc)) - irecp = system%irec + 1 + irecp = ireci + 1 do i = 1, plplenc_list%nenc - if (pl%levelg(plind1(i)) == irecp) pl%levelg(plind1(i)) = system%irec - if (pl%levelg(plind2(i)) == irecp) pl%levelg(plind2(i)) = system%irec + if (pl%levelg(plind1(i)) == irecp) pl%levelg(plind1(i)) = ireci + if (pl%levelg(plind2(i)) == irecp) pl%levelg(plind2(i)) = ireci end do do i = 1, pltpenc_list%nenc - if (pl%levelg(plind3(i)) == irecp) pl%levelg(plind3(i)) = system%irec - if (tp%levelg(tpind(i)) == irecp) tp%levelg(tpind(i)) = system%irec + if (pl%levelg(plind3(i)) == irecp) pl%levelg(plind3(i)) = ireci + if (tp%levelg(tpind(i)) == irecp) tp%levelg(tpind(i)) = ireci end do end associate - where(plplenc_list%level(1:plplenc_list%nenc) == irecp) plplenc_list%level(:) = system%irec - where(pltpenc_list%level(1:pltpenc_list%nenc) == irecp) pltpenc_list%level(:) = system%irec + if (plplenc_list%nenc > 0) where(plplenc_list%level(1:plplenc_list%nenc) == irecp) plplenc_list%level(:) = ireci + if (pltpenc_list%nenc > 0) where(pltpenc_list%level(1:pltpenc_list%nenc) == irecp) pltpenc_list%level(:) = ireci + system%irec = ireci end select end select end associate @@ -205,7 +207,7 @@ module recursive subroutine symba_step_recur_system(self, param, t, ireci) call pltpenc_list%collision_check(system, param, t+dtl, dtl, ireci) end if - call self%set_recur_levels() + call self%set_recur_levels(ireci) end do end select diff --git a/src/symba/symba_util.f90 b/src/symba/symba_util.f90 index 10fb36a2d..7a6f17cbf 100644 --- a/src/symba/symba_util.f90 +++ b/src/symba/symba_util.f90 @@ -59,13 +59,21 @@ module subroutine symba_util_resize_pltpenc(self, nrequested) ! Internals class(symba_pltpenc), allocatable :: enc_temp integer(I4B) :: nold + logical :: lmalloc - nold = size(self%status) + lmalloc = allocated(self%status) + if (lmalloc) then + nold = size(self%status) + else + nold = 0 + end if if (nrequested > nold) then - allocate(enc_temp, source=self) + if (lmalloc) allocate(enc_temp, source=self) call self%setup(2 * nrequested) - call self%copy(enc_temp) - deallocate(enc_temp) + if (lmalloc) then + call self%copy(enc_temp) + deallocate(enc_temp) + end if else self%status(nrequested+1:nold) = INACTIVE end if