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

Commit

Permalink
Browse files Browse the repository at this point in the history 
Changing back to ntp on outer loop and npl on inner loop
  • Loading branch information
@daminton
daminton committed Apr 14, 2021
1 parent 6e48ee5 commit d6c6369
Showing 1 changed file with 11 additions and 28 deletions.
39 changes: 11 additions & 28 deletions src/whm/whm_getacch.f90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -216,38 +216,21 @@ pure subroutine whm_getacch_ah3_tp(cb, pl, tp, xh)
integer(I4B) :: i, j
real(DP) :: rji2, irij3, fac
real(DP), dimension(NDIM) :: dx, acc
real(DP), dimension(:,:), allocatable :: aht

associate(ntp => tp%nbody, npl => pl%nbody, msun => cb%Gmass, GMpl => pl%Gmass, xht => tp%xh)

associate(ntp => tp%nbody, npl => pl%nbody, msun => cb%Gmass, GMpl => pl%Gmass, xht => tp%xh, aht => tp%ah)
if (ntp == 0) return
if (ntp > npl) then
allocate(aht, source=tp%ah)
do i = 1, ntp
acc(:) = 0.0_DP
!$omp simd private(dx,rji2,irij3,fac) reduction(-:acc)
do j = 1, npl
!$omp simd private(dx,rji2,irij3,fac) reduction(-:aht)
do i = 1, ntp
dx(:) = xht(:, i) - xh(:, j)
rji2 = dot_product(dx(:), dx(:))
irij3 = 1.0_DP / (rji2 * sqrt(rji2))
fac = GMpl(j) * irij3
aht(:, i) = aht(:, i) - fac * dx(:)
end do
end do
call move_alloc(aht, tp%ah)
else
do i = 1, ntp
acc(:) = 0.0_DP
!$omp simd private(dx,rji2,irij3,fac) reduction(-:acc)
do j = 1, npl
dx(:) = xht(:, i) - xh(:, j)
rji2 = dot_product(dx(:), dx(:))
irij3 = 1.0_DP / (rji2 * sqrt(rji2))
fac = GMpl(j) * irij3
acc(:) = acc(:) - fac * dx(:)
end do
tp%ah(:, i) = tp%ah(:, i) + acc(:)
dx(:) = xht(:, i) - xh(:, j)
rji2 = dot_product(dx(:), dx(:))
irij3 = 1.0_DP / (rji2 * sqrt(rji2))
fac = GMpl(j) * irij3
acc(:) = acc(:) - fac * dx(:)
end do
end if
aht(:, i) = aht(:, i) + acc(:)
end do
end associate
return
end subroutine whm_getacch_ah3_tp
Expand Down

0 comments on commit d6c6369

Please sign in to comment.