From 9bd9aca223399882f347ab4a858667d37ff8e1c6 Mon Sep 17 00:00:00 2001 From: David Minton Date: Fri, 19 May 2023 17:51:12 -0400 Subject: [PATCH] Set locality-spec for a bunch of do concurrents --- src/fraggle/fraggle_generate.f90 | 10 +++++++++- src/swiftest/swiftest_obl.f90 | 4 ++++ src/swiftest/swiftest_util.f90 | 23 +++++++++++++++++++++-- 3 files changed, 34 insertions(+), 3 deletions(-) diff --git a/src/fraggle/fraggle_generate.f90 b/src/fraggle/fraggle_generate.f90 index f1ba8d606..ee4466c11 100644 --- a/src/fraggle/fraggle_generate.f90 +++ b/src/fraggle/fraggle_generate.f90 @@ -299,7 +299,7 @@ module subroutine fraggle_generate_merge(self, nbody_system, param, t) real(DP), intent(in) :: t !! The time of the collision ! Internals - integer(I4B) :: i, j + integer(I4B) :: i real(DP), dimension(NDIM) :: L_spin_new, Ip, rot real(DP) :: rotmag, mass, volume, radius @@ -413,7 +413,11 @@ module subroutine fraggle_generate_pos_vec(collider, nbody_system, param, lfailu end do ! Randomly place the n>2 fragments inside their cloud until none are overlapping +#ifdef DOCONLOC + do concurrent(i = istart:nfrag, loverlap(i)) shared(loverlap, mass_rscale, u, phi, theta, lhitandrun) local(j, direction) +#else do concurrent(i = istart:nfrag, loverlap(i)) +#endif j = fragments%origin_body(i) ! Scale the cloud size @@ -452,7 +456,11 @@ module subroutine fraggle_generate_pos_vec(collider, nbody_system, param, lfailu ! Because body 1 and 2 are initialized near the original impactor positions, then if these bodies are still overlapping ! when the rest are not, we will randomly walk their position in space so as not to move them too far from their starting position if (all(.not.loverlap(istart:nfrag)) .and. any(loverlap(1:istart-1))) then +#ifdef DOCONLOC + do concurrent(i = 1:istart-1,loverlap(i)) shared(loverlap, u, theta, i) local(rwalk, dis) +#else do concurrent(i = 1:istart-1,loverlap(i)) +#endif dis = 0.1_DP * fragments%radius(i) * u(i)**(THIRD) rwalk(1) = fragments%rmag(i) * sin(theta(i)) * cos(phi(i)) rwalk(2) = fragments%rmag(i) * sin(theta(i)) * sin(phi(i)) diff --git a/src/swiftest/swiftest_obl.f90 b/src/swiftest/swiftest_obl.f90 index 7a6463677..85f88bc7d 100644 --- a/src/swiftest/swiftest_obl.f90 +++ b/src/swiftest/swiftest_obl.f90 @@ -147,7 +147,11 @@ module subroutine swiftest_obl_pot_system(self) associate(nbody_system => self, pl => self%pl, cb => self%cb) npl = self%pl%nbody if (.not. any(pl%lmask(1:npl))) return +#ifdef DOCONLOC + do concurrent (i = 1:npl, pl%lmask(i)) shared(oblpot_arr) +#else do concurrent (i = 1:npl, pl%lmask(i)) +#endif oblpot_arr(i) = swiftest_obl_pot_one(cb%Gmass, pl%Gmass(i), cb%j2rp2, cb%j4rp4, pl%rh(3,i), 1.0_DP / norm2(pl%rh(:,i))) end do nbody_system%oblpot = sum(oblpot_arr, pl%lmask(1:npl)) diff --git a/src/swiftest/swiftest_util.f90 b/src/swiftest/swiftest_util.f90 index 1784a1252..b3e43c7e2 100644 --- a/src/swiftest/swiftest_util.f90 +++ b/src/swiftest/swiftest_util.f90 @@ -1178,7 +1178,11 @@ module subroutine swiftest_util_get_energy_and_momentum_system(self, param) nbody_system%be_cb = -3*cb%Gmass * cb%mass / (5 * cb%radius) Lcborbit(:) = cb%mass * (cb%rb(:) .cross. cb%vb(:)) +#ifdef DOCONLOC + do concurrent (i = 1:npl, pl%lmask(i)) local(h) shared(Lplorbit, kepl) +#else do concurrent (i = 1:npl, pl%lmask(i)) +#endif h(:) = pl%rb(:,i) .cross. pl%vb(:,i) ! Angular momentum from orbit @@ -1194,7 +1198,11 @@ module subroutine swiftest_util_get_energy_and_momentum_system(self, param) ! For simplicity, we always assume that the rotation pole is the 3rd principal axis Lcbspin(:) = cb%Ip(3) * cb%mass * cb%radius**2 * cb%rot(:) +#ifdef DOCONLOC + do concurrent (i = 1:npl, pl%lmask(i)) shared(Lplspin, kespinpl) +#else do concurrent (i = 1:npl, pl%lmask(i)) +#endif ! Currently we assume that the rotation pole is the 3rd principal axis ! Angular momentum from spin Lplspin(:,i) = pl%mass(i) * pl%Ip(3,i) * pl%radius(i)**2 * pl%rot(:,i) @@ -1226,7 +1234,7 @@ module subroutine swiftest_util_get_energy_and_momentum_system(self, param) end if nbody_system%ke_orbit = 0.5_DP * (kecb + sum(kepl(1:npl), pl%lmask(1:npl))) - + do concurrent (j = 1:NDIM) nbody_system%L_orbit(j) = Lcborbit(j) + sum(Lplorbit(j,1:npl), pl%lmask(1:npl)) end do @@ -1271,7 +1279,11 @@ module subroutine swiftest_util_get_potential_energy_flat(npl, nplpl, k_plpl, lm pecb(1:npl) = 0.0_DP end where +#ifdef DOCONLOC + do concurrent(i = 1:npl, lmask(i)) shared(pecb) +#else do concurrent(i = 1:npl, lmask(i)) +#endif pecb(i) = -GMcb * mass(i) / norm2(rb(:,i)) end do @@ -1318,7 +1330,11 @@ module subroutine swiftest_util_get_potential_energy_triangular(npl, lmask, GMcb pecb(1:npl) = 0.0_DP end where +#ifdef DOCONLOC + do concurrent(i = 1:npl, lmask(i)) shared(pecb) +#else do concurrent(i = 1:npl, lmask(i)) +#endif pecb(i) = -GMcb * mass(i) / norm2(rb(:,i)) end do @@ -1329,7 +1345,11 @@ module subroutine swiftest_util_get_potential_energy_triangular(npl, lmask, GMcb !$omp reduction(+:pe) do i = 1, npl if (lmask(i)) then +#ifdef DOCONLOC + do concurrent(j = i+1:npl, lmask(i) .and. lmask(j)) shared(pepl) +#else do concurrent(j = i+1:npl, lmask(i) .and. lmask(j)) +#endif pepl(j) = - (Gmass(i) * mass(j)) / norm2(rb(:, i) - rb(:, j)) end do pe = pe + sum(pepl(i+1:npl), lmask(i+1:npl)) @@ -1523,7 +1543,6 @@ module subroutine swiftest_util_peri(n,m, r, v, atp, q, isperi) integer(I4B) :: i real(DP), dimension(n) :: e !! Temporary, just to make use of the xv2aeq subroutine real(DP) :: vdotr - character(len=NAMELEN) :: message do i = 1,n vdotr = dot_product(r(:,i),v(:,i))