diff --git a/Makefile.Defines b/Makefile.Defines index 70069bb71..291f2c604 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/examples/symba_swifter_comparison/1pl_1pl_encounter/cb.swiftest.in b/examples/symba_swifter_comparison/1pl_1pl_encounter/cb.swiftest.in index d0ae0ed15..4c5d87040 100644 Binary files a/examples/symba_swifter_comparison/1pl_1pl_encounter/cb.swiftest.in and b/examples/symba_swifter_comparison/1pl_1pl_encounter/cb.swiftest.in differ diff --git a/examples/symba_swifter_comparison/1pl_1pl_encounter/init_cond.py b/examples/symba_swifter_comparison/1pl_1pl_encounter/init_cond.py index f63022624..7600320c2 100755 --- a/examples/symba_swifter_comparison/1pl_1pl_encounter/init_cond.py +++ b/examples/symba_swifter_comparison/1pl_1pl_encounter/init_cond.py @@ -132,7 +132,7 @@ plfile = FortranFile(swiftest_pl, 'w') plfile.write_record(npl) -plfile.write_record(plid1) +plfile.write_record(np.array([plid1, plid2])) plfile.write_record(np.vstack([p_pl1[0],p_pl2[0]])) plfile.write_record(np.vstack([p_pl1[1],p_pl2[1]])) plfile.write_record(np.vstack([p_pl1[2],p_pl2[2]])) diff --git a/examples/symba_swifter_comparison/1pl_1pl_encounter/param.swiftest.in b/examples/symba_swifter_comparison/1pl_1pl_encounter/param.swiftest.in index a7f91ba33..1866557b2 100644 --- a/examples/symba_swifter_comparison/1pl_1pl_encounter/param.swiftest.in +++ b/examples/symba_swifter_comparison/1pl_1pl_encounter/param.swiftest.in @@ -5,7 +5,7 @@ DT 0.0006844626967830253 CB_IN cb.swiftest.in PL_IN pl.swiftest.in TP_IN tp.swiftest.in -IN_TYPE REAL8 +IN_TYPE ASCII ISTEP_OUT 1 ISTEP_DUMP 1 BIN_OUT bin.swiftest.dat diff --git a/examples/symba_swifter_comparison/1pl_1pl_encounter/pl.swifter.in b/examples/symba_swifter_comparison/1pl_1pl_encounter/pl.swifter.in index 9f0548fc1..0eb21018b 100644 --- a/examples/symba_swifter_comparison/1pl_1pl_encounter/pl.swifter.in +++ b/examples/symba_swifter_comparison/1pl_1pl_encounter/pl.swifter.in @@ -1,5 +1,5 @@ 3 ! Planet input file generated using init_cond.py -1 39.476926408897625196 +1 39.47692640889762629 0.0 0.0 0.0 0.0 0.0 0.0 2 0.00012002693582795244940133 0.010044724833237892 diff --git a/examples/symba_swifter_comparison/1pl_1pl_encounter/pl.swiftest.in b/examples/symba_swifter_comparison/1pl_1pl_encounter/pl.swiftest.in index 51f919531..19c6d6e3a 100644 Binary files a/examples/symba_swifter_comparison/1pl_1pl_encounter/pl.swiftest.in and b/examples/symba_swifter_comparison/1pl_1pl_encounter/pl.swiftest.in differ diff --git a/examples/symba_swifter_comparison/1pl_1pl_encounter/swiftest_vs_swifter.ipynb b/examples/symba_swifter_comparison/1pl_1pl_encounter/swiftest_vs_swifter.ipynb index dc1a9992f..57dd1934a 100644 --- a/examples/symba_swifter_comparison/1pl_1pl_encounter/swiftest_vs_swifter.ipynb +++ b/examples/symba_swifter_comparison/1pl_1pl_encounter/swiftest_vs_swifter.ipynb @@ -43,7 +43,7 @@ "output_type": "stream", "text": [ "Reading Swiftest file param.swiftest.in\n", - "Reading in time 6.845e-04\n", + "Reading in time 1.506e-01\n", "Creating Dataset\n" ] }, diff --git a/examples/symba_swifter_comparison/1pl_1pl_encounter/tp.swiftest.in b/examples/symba_swifter_comparison/1pl_1pl_encounter/tp.swiftest.in index 64bf92f74..573541ac9 100644 Binary files a/examples/symba_swifter_comparison/1pl_1pl_encounter/tp.swiftest.in and b/examples/symba_swifter_comparison/1pl_1pl_encounter/tp.swiftest.in differ diff --git a/python/swiftest/tests/bin2xr/swiftest/bin2xr_swiftest.py b/python/swiftest/tests/bin2xr/swiftest/bin2xr_swiftest.py index 9a377876a..0a430c654 100644 --- a/python/swiftest/tests/bin2xr/swiftest/bin2xr_swiftest.py +++ b/python/swiftest/tests/bin2xr/swiftest/bin2xr_swiftest.py @@ -1,5 +1,4 @@ import swiftest -import swiftest.io as swio -param_file = "param.swiftest.in" sim = swiftest.Simulation(param_file="param.swiftest.in") -ds = swio.swiftest2xr(sim.param) \ No newline at end of file +sim.bin2xr() +sim.ds \ No newline at end of file diff --git a/python/swiftest/tests/bin2xr/swiftest/param.swiftest.in b/python/swiftest/tests/bin2xr/swiftest/param.swiftest.in index 6a9e599f9..a7f91ba33 100644 --- a/python/swiftest/tests/bin2xr/swiftest/param.swiftest.in +++ b/python/swiftest/tests/bin2xr/swiftest/param.swiftest.in @@ -1,31 +1,31 @@ -! VERSION Swiftest parameter file converted from Swift -BIG_DISCARD NO -BIN_OUT bin.dat -CB_IN cb.swiftest.in -CHK_CLOSE YES -CHK_QMIN 0.00468 -CHK_QMIN_COORD HELIO -CHK_QMIN_RANGE 4.68e-03 100.0 -CHK_RMAX 100.0 -CHK_RMIN 0.00468 -DT 5.0 -DU2M 149597870700.0 -ENC_OUT enc.dat -ENERGY NO -EXTRA_FORCE NO -FRAGMENTATION NO -GR NO -ISTEP_DUMP 7305000 -ISTEP_OUT 7305000 -MU2KG 1.988409870698051e+30 -OUT_FORM EL -OUT_STAT UNKNOWN -PL_IN pl.swiftest.in -ROTATION NO -T0 0.0 -TIDES NO -TP_IN tp.swiftest.in -TSTOP 365250000000.0 -TU2S 86400 -YARKOVSKY NO -YORP NO +! Swiftest input file generated using init_cond.py +T0 0 +TSTOP 0.15 +DT 0.0006844626967830253 +CB_IN cb.swiftest.in +PL_IN pl.swiftest.in +TP_IN tp.swiftest.in +IN_TYPE REAL8 +ISTEP_OUT 1 +ISTEP_DUMP 1 +BIN_OUT bin.swiftest.dat +OUT_TYPE REAL8 +OUT_FORM XV +OUT_STAT REPLACE +CHK_CLOSE yes +CHK_RMIN 0.004650467260962157 +CHK_RMAX 1000.0 +CHK_EJECT 1000.0 +CHK_QMIN 0.004650467260962157 +CHK_QMIN_COORD HELIO +CHK_QMIN_RANGE 0.004650467260962157 1000.0 +ENC_OUT enc.swiftest.dat +EXTRA_FORCE no +BIG_DISCARD no +ROTATION no +GR no +MU2KG 1.988409870698051e+30 +DU2M 149597870700.0 +TU2S 31557600.0 +RHILL_PRESENT yes +MTINY 1e-12 diff --git a/python/swiftest/tests/bin2xr/swiftest/pl.swiftest.in b/python/swiftest/tests/bin2xr/swiftest/pl.swiftest.in index b7814ac98..d8da7a92a 100644 Binary files a/python/swiftest/tests/bin2xr/swiftest/pl.swiftest.in and b/python/swiftest/tests/bin2xr/swiftest/pl.swiftest.in differ diff --git a/python/swiftest/tests/bin2xr/swiftest/tp.swiftest.in b/python/swiftest/tests/bin2xr/swiftest/tp.swiftest.in index b3b8f423c..64bf92f74 100644 Binary files a/python/swiftest/tests/bin2xr/swiftest/tp.swiftest.in and b/python/swiftest/tests/bin2xr/swiftest/tp.swiftest.in differ diff --git a/src/symba/symba_encounter_check.f90 b/src/symba/symba_encounter_check.f90 index 1a9d8c68f..808ee2347 100644 --- a/src/symba/symba_encounter_check.f90 +++ b/src/symba/symba_encounter_check.f90 @@ -71,10 +71,11 @@ module function symba_encounter_check_pltpenc(self, system, dt, irec) result(lan 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 + integer(I4B) :: k real(DP), dimension(NDIM) :: xr, vr logical :: lencounter, isplpl real(DP) :: rlim2, rji2 + logical, dimension(:), allocatable :: lencmask lany_encounter = .false. if (self%nenc == 0) return @@ -90,40 +91,43 @@ module function symba_encounter_check_pltpenc(self, system, dt, irec) result(lan class is (symba_pl) select type(tp => system%tp) class is (symba_tp) - do concurrent(i = 1:self%nenc, self%status(i) == ACTIVE .and. self%level(i) == irec - 1) - associate(index_i => self%index1(i), index_j => self%index2(i)) + allocate(lencmask(self%nenc)) + lencmask(:) = (self%status(1:self%nenc) == ACTIVE) .and. (self%level(1:self%nenc) == irec - 1) + if (.not.any(lencmask(:))) return + associate(ind1 => self%index1, ind2 => self%index2) + do concurrent(k = 1:self%nenc, lencmask(k)) if (isplpl) then - xr(:) = pl%xh(:,index_j) - pl%xh(:,index_i) - vr(:) = pl%vb(:,index_j) - pl%vb(:,index_i) - call symba_encounter_check_one(xr(1), xr(2), xr(3), vr(1), vr(2), vr(3), pl%rhill(index_i), pl%rhill(index_j), dt, irec, lencounter, self%lvdotr(i)) + xr(:) = pl%xh(:,ind2(k)) - pl%xh(:,ind1(k)) + vr(:) = pl%vb(:,ind2(k)) - pl%vb(:,ind1(k)) + call symba_encounter_check_one(xr(1), xr(2), xr(3), vr(1), vr(2), vr(3), pl%rhill(ind1(k)), pl%rhill(ind2(k)), dt, irec, lencounter, self%lvdotr(k)) else - xr(:) = tp%xh(:,index_j) - pl%xh(:,index_i) - vr(:) = tp%vb(:,index_j) - pl%vb(:,index_i) - call symba_encounter_check_one(xr(1), xr(2), xr(3), vr(1), vr(2), vr(3), pl%rhill(index_i), 0.0_DP, dt, irec, lencounter, self%lvdotr(i)) + xr(:) = tp%xh(:,ind2(k)) - pl%xh(:,ind1(k)) + vr(:) = tp%vb(:,ind2(k)) - pl%vb(:,ind1(k)) + call symba_encounter_check_one(xr(1), xr(2), xr(3), vr(1), vr(2), vr(3), pl%rhill(ind1(k)), 0.0_DP, dt, irec, lencounter, self%lvdotr(k)) end if if (lencounter) then if (isplpl) then - rlim2 = (pl%radius(index_i) + pl%radius(index_j))**2 + rlim2 = (pl%radius(ind1(k)) + pl%radius(ind2(k)))**2 else - rlim2 = (pl%radius(index_i))**2 + rlim2 = (pl%radius(ind1(k)))**2 end if rji2 = dot_product(xr(:), xr(:))! Check to see if these are physically overlapping bodies first, which we should ignore if (rji2 > rlim2) then lany_encounter = .true. - pl%levelg(index_i) = irec - pl%levelm(index_i) = MAX(irec, pl%levelm(index_i)) + pl%levelg(ind1(k)) = irec + pl%levelm(ind1(k)) = MAX(irec, pl%levelm(ind1(k))) if (isplpl) then - pl%levelg(index_j) = irec - pl%levelm(index_j) = MAX(irec, pl%levelm(index_j)) + pl%levelg(ind2(k)) = irec + pl%levelm(ind2(k)) = MAX(irec, pl%levelm(ind2(k))) else - tp%levelg(index_j) = irec - tp%levelm(index_j) = MAX(irec, tp%levelm(index_j)) + tp%levelg(ind2(k)) = irec + tp%levelm(ind2(k)) = MAX(irec, tp%levelm(ind2(k))) end if - self%level(i) = irec + self%level(k) = irec end if end if - end associate - end do + end do + end associate end select end select