Skip to content

Replace OMP reductions #137

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
marvinfriede merged 3 commits into from
Apr 8, 2025
Merged

Replace OMP reductions #137

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
3 commits
Select commit Hold shift + click to select a range
06da3a1
Remove unused vars and ignore wrap files
marvinfriede Apr 8, 2025
57c807f
Replace OMP reductions
marvinfriede Apr 8, 2025
5f7f5ad
Merge branch 'main' into reductions
marvinfriede Apr 8, 2025
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
3 changes: 1 addition & 2 deletions app/driver.f90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -71,8 +71,7 @@ subroutine run_driver(config, error)
real(wp), allocatable :: s9
real(wp) :: energy
character(len=:), allocatable :: output
integer :: stat, unit, idx
logical :: exist
integer :: unit
type(citation_type) :: citation, param_citation

if (config%verbosity > 1) then
Expand Down
115 changes: 81 additions & 34 deletions src/dftd3/damping/atm.f90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -119,13 +119,21 @@ subroutine get_atm_dispersion_energy(mol, trans, cutoff, s9, rs9, alp, rvdw, c6,
real(wp) :: r0ij, r0jk, r0ik, r0, r1, r2, r3, r5, rr, fdmp, ang
real(wp) :: cutoff2, c9, dE

! Thread-private arrays for reduction
! Set to 0 explicitly as the shared variants are potentially non-zero (inout)
real(wp), allocatable :: energy_local(:)

cutoff2 = cutoff*cutoff

!$omp parallel do schedule(runtime) default(none) reduction(+:energy) &
!$omp parallel default(none) &
!$omp shared(mol, trans, c6, s9, rs9, alp, rvdw, cutoff2) &
!$omp private(iat, jat, kat, izp, jzp, kzp, jtr, ktr, vij, vjk, vik, &
!$omp& r2ij, r2jk, r2ik, c6ij, c6jk, c6ik, triple, r0ij, r0jk, r0ik, r0, &
!$omp& r1, r2, r3, r5, rr, fdmp, ang, c9, dE)
!$omp& r1, r2, r3, r5, rr, fdmp, ang, c9, dE) &
!$omp shared(energy) &
!$omp private(energy_local)
allocate(energy_local(size(energy, 1)), source=0.0_wp)
!$omp do schedule(runtime)
do iat = 1, mol%nat
izp = mol%id(iat)
do jat = 1, iat
Expand Down Expand Up @@ -164,15 +172,21 @@ subroutine get_atm_dispersion_energy(mol, trans, cutoff, s9, rs9, alp, rvdw, c6,

rr = ang*fdmp

dE = rr * c9 * triple
energy(iat) = energy(iat) - dE/3
energy(jat) = energy(jat) - dE/3
energy(kat) = energy(kat) - dE/3
dE = rr * c9 * triple / 3.0_wp
energy_local(iat) = energy_local(iat) - dE
energy_local(jat) = energy_local(jat) - dE
energy_local(kat) = energy_local(kat) - dE
end do
end do
end do
end do
end do
!$omp end do
!$omp critical (get_atm_dispersion_energy_)
energy(:) = energy(:) + energy_local(:)
!$omp end critical (get_atm_dispersion_energy_)
deallocate(energy_local)
!$omp end parallel

end subroutine get_atm_dispersion_energy

Expand Down Expand Up @@ -225,15 +239,28 @@ subroutine get_atm_dispersion_derivs(mol, trans, cutoff, s9, rs9, alp, rvdw, c6,
real(wp) :: r0ij, r0jk, r0ik, r0, r1, r2, r3, r5, rr, fdmp, dfdmp, ang, dang
real(wp) :: cutoff2, c9, dE, dGij(3), dGjk(3), dGik(3), dS(3, 3)

! Thread-private arrays for reduction
! Set to 0 explicitly as the shared variants are potentially non-zero (inout)
real(wp), allocatable :: energy_local(:)
real(wp), allocatable :: dEdcn_local(:)
real(wp), allocatable :: gradient_local(:, :)
real(wp), allocatable :: sigma_local(:, :)

cutoff2 = cutoff*cutoff

!$omp parallel do schedule(runtime) default(none) &
!$omp reduction(+:energy, gradient, sigma, dEdcn) &
!$omp parallel default(none) &
!$omp shared(mol, trans, c6, s9, rs9, alp, rvdw, cutoff2, dc6dcn) &
!$omp private(iat, jat, kat, izp, jzp, kzp, jtr, ktr, vij, vjk, vik, &
!$omp& r2ij, r2jk, r2ik, c6ij, c6jk, c6ik, triple, r0ij, r0jk, r0ik, r0, &
!$omp& r1, r2, r3, r5, rr, fdmp, dfdmp, ang, dang, c9, dE, dGij, dGjk, &
!$omp& dGik, dS)
!$omp& dGik, dS) &
!$omp shared(energy, gradient, sigma, dEdcn) &
!$omp private(energy_local, gradient_local, sigma_local, dEdcn_local)
allocate(energy_local(size(energy, 1)), source=0.0_wp)
allocate(dEdcn_local(size(dEdcn, 1)), source=0.0_wp)
allocate(gradient_local(size(gradient, 1), size(gradient, 2)), source=0.0_wp)
allocate(sigma_local(size(sigma, 1), size(sigma, 2)), source=0.0_wp)
!$omp do schedule(runtime)
do iat = 1, mol%nat
izp = mol%id(iat)
do jat = 1, iat
Expand Down Expand Up @@ -296,31 +323,43 @@ subroutine get_atm_dispersion_derivs(mol, trans, cutoff, s9, rs9, alp, rvdw, c6,
dGjk(:) = c9 * (-dang * fdmp + ang * dfdmp) / r2jk * vjk

dE = rr * c9 * triple
energy(iat) = energy(iat) - dE/3
energy(jat) = energy(jat) - dE/3
energy(kat) = energy(kat) - dE/3
energy_local(iat) = energy_local(iat) - dE/3.0_wp
energy_local(jat) = energy_local(jat) - dE/3.0_wp
energy_local(kat) = energy_local(kat) - dE/3.0_wp

gradient(:, iat) = gradient(:, iat) - dGij - dGik
gradient(:, jat) = gradient(:, jat) + dGij - dGjk
gradient(:, kat) = gradient(:, kat) + dGik + dGjk
gradient_local(:, iat) = gradient_local(:, iat) - dGij - dGik
gradient_local(:, jat) = gradient_local(:, jat) + dGij - dGjk
gradient_local(:, kat) = gradient_local(:, kat) + dGik + dGjk

dS(:, :) = spread(dGij, 1, 3) * spread(vij, 2, 3)&
& + spread(dGik, 1, 3) * spread(vik, 2, 3)&
& + spread(dGjk, 1, 3) * spread(vjk, 2, 3)

sigma(:, :) = sigma + dS * triple
sigma_local(:, :) = sigma_local + dS * triple

dEdcn(iat) = dEdcn(iat) - dE * 0.5_wp &
dEdcn_local(iat) = dEdcn_local(iat) - dE * 0.5_wp &
& * (dc6dcn(iat, jat) / c6ij + dc6dcn(iat, kat) / c6ik)
dEdcn(jat) = dEdcn(jat) - dE * 0.5_wp &
dEdcn_local(jat) = dEdcn_local(jat) - dE * 0.5_wp &
& * (dc6dcn(jat, iat) / c6ij + dc6dcn(jat, kat) / c6jk)
dEdcn(kat) = dEdcn(kat) - dE * 0.5_wp &
dEdcn_local(kat) = dEdcn_local(kat) - dE * 0.5_wp &
& * (dc6dcn(kat, iat) / c6ik + dc6dcn(kat, jat) / c6jk)
end do
end do
end do
end do
end do
!$omp end do
!$omp critical (get_atm_dispersion_derivs_)
energy(:) = energy(:) + energy_local(:)
dEdcn(:) = dEdcn(:) + dEdcn_local(:)
gradient(:, :) = gradient(:, :) + gradient_local(:, :)
sigma(:, :) = sigma(:, :) + sigma_local(:, :)
!$omp end critical (get_atm_dispersion_derivs_)
deallocate(energy_local)
deallocate(dEdcn_local)
deallocate(gradient_local)
deallocate(sigma_local)
!$omp end parallel

end subroutine get_atm_dispersion_derivs

Expand Down Expand Up @@ -361,14 +400,22 @@ subroutine get_atm_pairwise_dispersion(mol, trans, cutoff, s9, rs9, alp, rvdw, c
real(wp) :: r0ij, r0jk, r0ik, r0, r1, r2, r3, r5, rr, fdmp, ang
real(wp) :: cutoff2, c9, dE

! Thread-private arrays for reduction
! Set to 0 explicitly as the shared variants are potentially non-zero (inout)
real(wp), allocatable :: energy_local(:, :)

if (abs(s9) < epsilon(1.0_wp)) return
cutoff2 = cutoff*cutoff

!$omp parallel do schedule(runtime) default(none) shared(energy) &
!$omp parallel default(none) &
!$omp shared(mol, trans, c6, cutoff2, s9, rs9, alp, rvdw) &
!$omp private(iat, jat, kat, izp, jzp, kzp, jtr, ktr, vij, vjk, vik, &
!$omp& r2ij, r2jk, r2ik, c6ij, c6jk, c6ik, triple, r0ij, r0jk, r0ik, r0, &
!$omp& r1, r2, r3, r5, rr, fdmp, ang, c9, dE)
!$omp& r1, r2, r3, r5, rr, fdmp, ang, c9, dE) &
!$omp shared(energy) &
!$omp private(energy_local)
allocate(energy_local(size(energy, 1), size(energy, 2)), source=0.0_wp)
!$omp do schedule(runtime)
do iat = 1, mol%nat
izp = mol%id(iat)
do jat = 1, iat
Expand Down Expand Up @@ -407,24 +454,24 @@ subroutine get_atm_pairwise_dispersion(mol, trans, cutoff, s9, rs9, alp, rvdw, c

rr = ang*fdmp

dE = rr * c9 * triple/6
!$omp atomic
energy(jat, iat) = energy(jat, iat) - dE
!$omp atomic
energy(kat, iat) = energy(kat, iat) - dE
!$omp atomic
energy(iat, jat) = energy(iat, jat) - dE
!$omp atomic
energy(kat, jat) = energy(kat, jat) - dE
!$omp atomic
energy(iat, kat) = energy(iat, kat) - dE
!$omp atomic
energy(jat, kat) = energy(jat, kat) - dE
dE = rr * c9 * triple / 6.0_wp
energy_local(jat, iat) = energy_local(jat, iat) - dE
energy_local(kat, iat) = energy_local(kat, iat) - dE
energy_local(iat, jat) = energy_local(iat, jat) - dE
energy_local(kat, jat) = energy_local(kat, jat) - dE
energy_local(iat, kat) = energy_local(iat, kat) - dE
energy_local(jat, kat) = energy_local(jat, kat) - dE
end do
end do
end do
end do
end do
!$omp end do
!$omp critical (get_atm_pairwise_dispersion_)
energy(:, :) = energy(:, :) + energy_local(:, :)
!$omp end critical (get_atm_pairwise_dispersion_)
deallocate(energy_local)
!$omp end parallel

end subroutine get_atm_pairwise_dispersion

Expand Down
93 changes: 72 additions & 21 deletions src/dftd3/damping/mzero.f90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -164,14 +164,22 @@ subroutine get_dispersion_energy(self, mol, trans, cutoff, rvdw, r4r2, c6, energ
real(wp) :: vec(3), r2, r1, r6, r8, t6, t8, f6, f8, alp6, alp8
real(wp) :: edisp, cutoff2, r0ij, rrij, c6ij, dE

! Thread-private array for reduction
! Set to 0 explicitly as the shared variants are potentially non-zero (inout)
real(wp), allocatable :: energy_local(:)

cutoff2 = cutoff*cutoff
alp6 = self%alp
alp8 = self%alp + 2.0_wp

!$omp parallel do schedule(runtime) default(none) reduction(+:energy) &
!$omp parallel default(none) &
!$omp shared(mol, self, c6, trans, cutoff2, alp6, alp8, rvdw, r4r2) &
!$omp private(iat, jat, izp, jzp, jtr, vec, r2, r1, r6, r8, t6, t8, f6, &
!$omp& f8, edisp, r0ij, rrij, c6ij, dE)
!$omp& f8, edisp, r0ij, rrij, c6ij, dE) &
!$omp shared(energy) &
!$omp private(energy_local)
allocate(energy_local(size(energy, 1)), source=0.0_wp)
!$omp do schedule(runtime)
do iat = 1, mol%nat
izp = mol%id(iat)
do jat = 1, iat
Expand All @@ -198,13 +206,19 @@ subroutine get_dispersion_energy(self, mol, trans, cutoff, rvdw, r4r2, c6, energ

dE = -c6ij*edisp * 0.5_wp

energy(iat) = energy(iat) + dE
energy_local(iat) = energy_local(iat) + dE
if (iat /= jat) then
energy(jat) = energy(jat) + dE
energy_local(jat) = energy_local(jat) + dE
end if
end do
end do
end do
!$omp end do
!$omp critical (get_dispersion_energy_)
energy(:) = energy(:) + energy_local(:)
!$omp end critical (get_dispersion_energy_)
deallocate(energy_local)
!$omp end parallel

end subroutine get_dispersion_energy

Expand Down Expand Up @@ -253,15 +267,28 @@ subroutine get_dispersion_derivs(self, mol, trans, cutoff, rvdw, r4r2, c6, dc6dc
real(wp) :: vec(3), r2, r1, r6, r8, t6, t8, d6, d8, f6, f8, alp6, alp8
real(wp) :: edisp, gdisp, cutoff2, r0ij, rrij, c6ij, dE, dG(3), dS(3, 3)

! Thread-private arrays for reduction
! Set to 0 explicitly as the shared variants are potentially non-zero (inout)
real(wp), allocatable :: energy_local(:)
real(wp), allocatable :: dEdcn_local(:)
real(wp), allocatable :: gradient_local(:, :)
real(wp), allocatable :: sigma_local(:, :)

cutoff2 = cutoff*cutoff
alp6 = self%alp
alp8 = self%alp + 2.0_wp

!$omp parallel do schedule(runtime) default(none) &
!$omp reduction(+:energy, gradient, sigma, dEdcn) &
!$omp parallel default(none) &
!$omp shared(mol, self, c6, dc6dcn, trans, cutoff2, alp6, alp8, rvdw, r4r2) &
!$omp private(iat, jat, izp, jzp, jtr, vec, r2, r1, r6, r8, t6, t8, d6, &
!$omp& d8, f6, f8, edisp, gdisp, r0ij, rrij, c6ij, dE, dG, dS)
!$omp& d8, f6, f8, edisp, gdisp, r0ij, rrij, c6ij, dE, dG, dS) &
!$omp shared(energy, gradient, sigma, dEdcn) &
!$omp private(energy_local, gradient_local, sigma_local, dEdcn_local)
allocate(energy_local(size(energy, 1)), source=0.0_wp)
allocate(dEdcn_local(size(dEdcn, 1)), source=0.0_wp)
allocate(gradient_local(size(gradient, 1), size(gradient, 2)), source=0.0_wp)
allocate(sigma_local(size(sigma, 1), size(sigma, 2)), source=0.0_wp)
!$omp do schedule(runtime)
do iat = 1, mol%nat
izp = mol%id(iat)
do jat = 1, iat
Expand Down Expand Up @@ -296,19 +323,31 @@ subroutine get_dispersion_derivs(self, mol, trans, cutoff, rvdw, r4r2, c6, dc6dc
dG(:) = -c6ij*gdisp*vec
dS(:, :) = spread(dG, 1, 3) * spread(vec, 2, 3) * 0.5_wp

energy(iat) = energy(iat) + dE
dEdcn(iat) = dEdcn(iat) - dc6dcn(iat, jat) * edisp
sigma(:, :) = sigma + dS
energy_local(iat) = energy_local(iat) + dE
dEdcn_local(iat) = dEdcn_local(iat) - dc6dcn(iat, jat) * edisp
sigma_local(:, :) = sigma_local + dS
if (iat /= jat) then
energy(jat) = energy(jat) + dE
dEdcn(jat) = dEdcn(jat) - dc6dcn(jat, iat) * edisp
gradient(:, iat) = gradient(:, iat) + dG
gradient(:, jat) = gradient(:, jat) - dG
sigma(:, :) = sigma + dS
energy_local(jat) = energy_local(jat) + dE
dEdcn_local(jat) = dEdcn_local(jat) - dc6dcn(jat, iat) * edisp
gradient_local(:, iat) = gradient_local(:, iat) + dG
gradient_local(:, jat) = gradient_local(:, jat) - dG
sigma_local(:, :) = sigma_local + dS
end if
end do
end do
end do
!$omp end do
!$omp critical (get_dispersion_derivs_)
energy(:) = energy(:) + energy_local(:)
dEdcn(:) = dEdcn(:) + dEdcn_local(:)
gradient(:, :) = gradient(:, :) + gradient_local(:, :)
sigma(:, :) = sigma(:, :) + sigma_local(:, :)
!$omp end critical (get_dispersion_derivs_)
deallocate(energy_local)
deallocate(dEdcn_local)
deallocate(gradient_local)
deallocate(sigma_local)
!$omp end parallel

end subroutine get_dispersion_derivs

Expand Down Expand Up @@ -390,14 +429,22 @@ subroutine get_pairwise_dispersion2(self, mol, trans, cutoff, rvdw, r4r2, c6, en
real(wp) :: vec(3), r2, r1, r6, r8, t6, t8, f6, f8, alp6, alp8
real(wp) :: edisp, cutoff2, r0ij, rrij, c6ij, dE

! Thread-private array for reduction
! Set to 0 explicitly as the shared variants are potentially non-zero (inout)
real(wp), allocatable :: energy_local(:, :)

cutoff2 = cutoff*cutoff
alp6 = self%alp
alp8 = self%alp + 2.0_wp

!$omp parallel do schedule(runtime) default(none) shared(energy) &
!$omp parallel default(none) &
!$omp shared(mol, self, c6, trans, cutoff2, alp6, alp8, rvdw, r4r2) &
!$omp private(iat, jat, izp, jzp, jtr, vec, r2, r1, r6, r8, t6, t8, f6, &
!$omp& f8, edisp, r0ij, rrij, c6ij, dE)
!$omp& f8, edisp, r0ij, rrij, c6ij, dE) &
!$omp shared(energy) &
!$omp private(energy_local)
allocate(energy_local(size(energy, 1), size(energy, 2)), source=0.0_wp)
!$omp do schedule(runtime)
do iat = 1, mol%nat
izp = mol%id(iat)
do jat = 1, iat
Expand All @@ -424,15 +471,19 @@ subroutine get_pairwise_dispersion2(self, mol, trans, cutoff, rvdw, r4r2, c6, en

dE = -c6ij*edisp * 0.5_wp

!$omp atomic
energy(jat, iat) = energy(jat, iat) + dE
energy_local(jat, iat) = energy_local(jat, iat) + dE
if (iat /= jat) then
!$omp atomic
energy(iat, jat) = energy(iat, jat) + dE
energy_local(iat, jat) = energy_local(iat, jat) + dE
end if
end do
end do
end do
!$omp end do
!$omp critical (get_pairwise_dispersion2_)
energy(:, :) = energy(:, :) + energy_local(:, :)
!$omp end critical (get_pairwise_dispersion2_)
deallocate(energy_local)
!$omp end parallel

end subroutine get_pairwise_dispersion2

Expand Down
Loading
Loading