Skip to content
This repository was archived by the owner on Aug 28, 2024. It is now read-only.

Commit

Permalink
Set locality-spec for a bunch of do concurrents
Browse files Browse the repository at this point in the history
  • Loading branch information
daminton committed May 19, 2023
1 parent 7ceea8c commit 9bd9aca
Show file tree
Hide file tree
Showing 3 changed files with 34 additions and 3 deletions.
10 changes: 9 additions & 1 deletion src/fraggle/fraggle_generate.f90
Original file line number Diff line number Diff line change
Expand Up @@ -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

Expand Down Expand Up @@ -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
Expand Down Expand Up @@ -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))
Expand Down
4 changes: 4 additions & 0 deletions src/swiftest/swiftest_obl.f90
Original file line number Diff line number Diff line change
Expand Up @@ -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))
Expand Down
23 changes: 21 additions & 2 deletions src/swiftest/swiftest_util.f90
Original file line number Diff line number Diff line change
Expand Up @@ -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
Expand All @@ -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)
Expand Down Expand Up @@ -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
Expand Down Expand Up @@ -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

Expand Down Expand Up @@ -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

Expand All @@ -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))
Expand Down Expand Up @@ -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))
Expand Down

0 comments on commit 9bd9aca

Please sign in to comment.