From 723dc4ef407be0ff189e145181f9ac1f48248a35 Mon Sep 17 00:00:00 2001 From: Paolo Giannozzi Date: Tue, 16 Jan 2018 22:13:52 +0100 Subject: [PATCH] Routine "gradrho" moved to gradutils, with new name "fft_gradient_g2r", while previous "fft_gradient" becomes "fft_gradient_r2r". Routine "grad_dot" moved to gradutils, with new name "fft_graddot" and removal of useless variable in the list of argument. --- LR_Modules/compute_vsgga.f90 | 3 +- LR_Modules/setup_dgc.f90 | 5 +- Modules/gradutils.f90 | 140 ++++++++++++++++++++++++++++++++++- PHonon/Gamma/cg_setupdgc.f90 | 2 +- PP/src/elf.f90 | 2 +- PW/src/gradcorr.f90 | 131 +------------------------------- PW/src/loc_scdm.f90 | 2 +- PW/src/stres_gradcorr.f90 | 2 +- PW/src/v_of_rho.f90 | 6 +- 9 files changed, 150 insertions(+), 143 deletions(-) diff --git a/LR_Modules/compute_vsgga.f90 b/LR_Modules/compute_vsgga.f90 index ad12f2b56..8e134c636 100644 --- a/LR_Modules/compute_vsgga.f90 +++ b/LR_Modules/compute_vsgga.f90 @@ -13,7 +13,6 @@ SUBROUTINE compute_vsgga( rhoout, grho, vsgga ) USE constants, ONLY : e2 USE kinds, ONLY : DP USE gvect, ONLY : ngm, g - USE cell_base, ONLY : alat USE noncollin_module, ONLY : noncolin, nspin_gga USE funct, ONLY : gcxc, gcx_spin, gcc_spin, & gcc_spin_more, dft_is_gradient, get_igcc @@ -154,7 +153,7 @@ SUBROUTINE compute_vsgga( rhoout, grho, vsgga ) ! DO is = 1, nspin_gga ! - CALL grad_dot( dfftp, h(1,1,is), g, alat, dh ) + CALL fft_graddot( dfftp, h(1,1,is), g, dh ) ! vaux(:,is) = vaux(:,is) - dh(:) ! diff --git a/LR_Modules/setup_dgc.f90 b/LR_Modules/setup_dgc.f90 index 258553072..a2a3cdb27 100644 --- a/LR_Modules/setup_dgc.f90 +++ b/LR_Modules/setup_dgc.f90 @@ -79,8 +79,7 @@ subroutine setup_dgc ! rhogout(:,is) = psic(dfftp%nl(:)) ! - ! - CALL gradrho(dfftp, rhogout(1,is), g, grho(1,1,is) ) + CALL fft_gradient_g2r(dfftp, rhogout(1,is), g, grho(1,1,is) ) ! END DO DEALLOCATE(rhogout) @@ -95,7 +94,7 @@ subroutine setup_dgc enddo endif do is = 1, nspin_gga - call gradrho (dfftp, rho%of_g (1, is), g, grho (1, 1, is) ) + call fft_gradient_g2r (dfftp, rho%of_g (1, is), g, grho (1, 1, is) ) enddo END IF diff --git a/Modules/gradutils.f90 b/Modules/gradutils.f90 index fdb16e5d0..419e4282b 100644 --- a/Modules/gradutils.f90 +++ b/Modules/gradutils.f90 @@ -26,13 +26,13 @@ SUBROUTINE external_gradient( a, grada ) REAL( DP ), INTENT(OUT) :: grada( 3, dfftp%nnr ) ! A in real space, grad(A) in real space - CALL ffT_gradient( dfftp, a, g, grada ) + CALL fft_gradient_r2r( dfftp, a, g, grada ) RETURN END SUBROUTINE external_gradient !---------------------------------------------------------------------------- -SUBROUTINE fft_gradient( dfft, a, g, ga ) +SUBROUTINE fft_gradient_r2r( dfft, a, g, ga ) !---------------------------------------------------------------------------- ! ! ... Calculates ga = \grad a @@ -96,8 +96,142 @@ SUBROUTINE fft_gradient( dfft, a, g, ga ) ! RETURN ! -END SUBROUTINE fft_gradient +END SUBROUTINE fft_gradient_r2r !-------------------------------------------------------------------- +! +!---------------------------------------------------------------------------- +SUBROUTINE fft_gradient_g2r( dfft, a, g, ga ) + !---------------------------------------------------------------------------- + ! + ! ... Calculates ga = \grad a - like fft_gradient with a(G) instead of a(r) + ! ... input : dfft FFT descriptor + ! ... a(:) a(G), a complex function in G-space + ! ... g(3,:) G-vectors, in 2pi/a units + ! ... output: ga(3,:) \grad a, real, on the real-space FFT grid + ! + USE cell_base, ONLY : tpiba + USE kinds, ONLY : DP + USE control_flags, ONLY : gamma_only + USE fft_interfaces,ONLY : invfft + USE fft_types, ONLY : fft_type_descriptor + ! + IMPLICIT NONE + ! + TYPE(fft_type_descriptor),INTENT(IN) :: dfft + COMPLEX(DP), INTENT(IN) :: a(dfft%ngm) + REAL(DP), INTENT(IN) :: g(3,dfft%ngm) + REAL(DP), INTENT(OUT) :: ga(3,dfft%nnr) + ! + INTEGER :: ipol + COMPLEX(DP), ALLOCATABLE :: gaux(:) + ! + ! + ALLOCATE( gaux( dfft%nnr ) ) + ! + ! ... multiply by (iG) to get (\grad_ipol a)(G) ... + ! + ga(:,:) = 0.D0 + ! + DO ipol = 1, 3 + ! + gaux(:) = (0.0_dp,0.0_dp) + ! + gaux(dfft%nl(:)) = g(ipol,:) * CMPLX( -AIMAG(a(:)), REAL(a(:)), kind=DP) + ! + IF ( gamma_only ) THEN + ! + gaux(dfft%nlm(:)) = CMPLX( REAL( gaux(dfft%nl(:)) ), & + -AIMAG( gaux(dfft%nl(:)) ), kind=DP) + ! + END IF + ! + ! ... bring back to R-space, (\grad_ipol a)(r) ... + ! + CALL invfft ('Rho', gaux, dfft) + ! + ! ...and add the factor 2\pi/a missing in the definition of G + ! + ga(ipol,:) = ga(ipol,:) + tpiba * REAL( gaux(:) ) + ! + END DO + ! + DEALLOCATE( gaux ) + ! + RETURN + ! +END SUBROUTINE fft_gradient_g2r + +!---------------------------------------------------------------------------- +SUBROUTINE fft_graddot( dfft, a, g, da ) + !---------------------------------------------------------------------------- + ! + ! ... Calculates da = \sum_i \grad_i a_i in R-space + ! ... input : dfft FFT descriptor + ! ... a(3,:) a real function on the real-space FFT grid + ! ... g(3,:) G-vectors, in 2pi/a units + ! ... output: ga(:) \sum_i \grad_i a_i, real, on the real-space FFT grid + ! + USE cell_base, ONLY : tpiba + USE kinds, ONLY : DP + USE control_flags, ONLY : gamma_only + USE fft_interfaces,ONLY : fwfft, invfft + USE fft_types, ONLY : fft_type_descriptor + ! + IMPLICIT NONE + ! + TYPE(fft_type_descriptor),INTENT(IN) :: dfft + REAL(DP), INTENT(IN) :: a(3,dfft%nnr), g(3,dfft%ngm) + REAL(DP), INTENT(OUT) :: da(dfft%nnr) + ! + INTEGER :: n, ipol + COMPLEX(DP), ALLOCATABLE :: aux(:), gaux(:) + ! + ! + ALLOCATE( aux(dfft%nnr), gaux(dfft%nnr) ) + ! + gaux(:) = (0.0_dp,0.0_dp) + ! + DO ipol = 1, 3 + ! + aux = CMPLX( a(ipol,:), 0.0_dp, kind=DP) + ! + ! ... bring a(ipol,r) to G-space, a(G) ... + ! + CALL fwfft ('Rho', aux, dfft) + ! + DO n = 1, dfft%ngm + ! + gaux(dfft%nl(n)) = gaux(dfft%nl(n)) + g(ipol,n) * & + CMPLX( -AIMAG( aux(dfft%nl(n)) ), & + REAL( aux(dfft%nl(n)) ), kind=DP) + ! + END DO + ! + END DO + ! + IF ( gamma_only ) THEN + ! + DO n = 1, dfft%ngm + ! + gaux(dfft%nlm(n)) = CONJG( gaux(dfft%nl(n)) ) + ! + END DO + ! + END IF + ! + ! ... bring back to R-space, (\grad_ipol a)(r) ... + ! + CALL invfft ('Rho', gaux, dfft) + ! + ! ... add the factor 2\pi/a missing in the definition of G and sum + ! + da(:) = tpiba * REAL( gaux(:) ) + ! + DEALLOCATE( aux, gaux ) + ! + RETURN + ! +END SUBROUTINE fft_graddot !-------------------------------------------------------------------- ! Routines computing laplacian via FFT diff --git a/PHonon/Gamma/cg_setupdgc.f90 b/PHonon/Gamma/cg_setupdgc.f90 index d8865f5c9..9c66922c4 100644 --- a/PHonon/Gamma/cg_setupdgc.f90 +++ b/PHonon/Gamma/cg_setupdgc.f90 @@ -54,7 +54,7 @@ SUBROUTINE cg_setupdgc ENDDO ENDIF DO is=1,nspin - CALL gradrho (dfftp, rho%of_g(1,is), g, grho(1,1,is)) + CALL fft_gradient_g2r (dfftp, rho%of_g(1,is), g, grho(1,1,is)) ENDDO ! IF (nspin==1) THEN diff --git a/PP/src/elf.f90 b/PP/src/elf.f90 index 49a0afa5c..3577f96c7 100644 --- a/PP/src/elf.f90 +++ b/PP/src/elf.f90 @@ -210,7 +210,7 @@ SUBROUTINE do_rdg (rdg) ENDDO ! gradient of rho - CALL gradrho(dfftp, rho%of_g(1,1), g, grho) + CALL fft_gradient_g2r(dfftp, rho%of_g(1,1), g, grho) ! calculate rdg DO i = 1, dfftp%nnr diff --git a/PW/src/gradcorr.f90 b/PW/src/gradcorr.f90 index 9af1fd4ee..c8e4e7f17 100644 --- a/PW/src/gradcorr.f90 +++ b/PW/src/gradcorr.f90 @@ -14,7 +14,7 @@ SUBROUTINE gradcorr( rho, rhog, rho_core, rhog_core, etxc, vtxc, v ) USE kinds, ONLY : DP USE gvect, ONLY : ngm, g USE lsda_mod, ONLY : nspin - USE cell_base, ONLY : omega, alat + USE cell_base, ONLY : omega USE funct, ONLY : gcxc, gcx_spin, gcc_spin, igcc_is_lyp, & gcc_spin_more, dft_is_gradient, get_igcc USE spin_orb, ONLY : domag @@ -99,7 +99,7 @@ SUBROUTINE gradcorr( rho, rhog, rho_core, rhog_core, etxc, vtxc, v ) rhoout(:,is) = fac * rho_core(:) + rhoout(:,is) rhogsum(:,is) = fac * rhog_core(:) + rhogsum(:,is) ! - CALL gradrho( dfftp, rhogsum(1,is), g, grho(1,1,is) ) + CALL fft_gradient_g2r( dfftp, rhogsum(1,is), g, grho(1,1,is) ) ! END DO ! @@ -253,7 +253,7 @@ SUBROUTINE gradcorr( rho, rhog, rho_core, rhog_core, etxc, vtxc, v ) ! DO is = 1, nspin0 ! - CALL grad_dot( dfftp, h(1,1,is), g, alat, dh ) + CALL fft_graddot( dfftp, h(1,1,is), g, dh ) ! v(:,is) = v(:,is) - dh(:) ! @@ -292,128 +292,3 @@ SUBROUTINE gradcorr( rho, rhog, rho_core, rhog_core, etxc, vtxc, v ) RETURN ! END SUBROUTINE gradcorr -! -!---------------------------------------------------------------------------- -SUBROUTINE gradrho( dfft, a, g, ga ) - !---------------------------------------------------------------------------- - ! - ! ... Calculates ga = \grad a in R-space (a is in G-space) - ! - USE cell_base, ONLY : tpiba - USE kinds, ONLY : DP - USE control_flags, ONLY : gamma_only - USE fft_interfaces,ONLY : invfft - USE fft_types, ONLY : fft_type_descriptor - ! - IMPLICIT NONE - ! - TYPE(fft_type_descriptor),INTENT(IN) :: dfft - COMPLEX(DP), INTENT(IN) :: a(dfft%ngm) - REAL(DP), INTENT(IN) :: g(3,dfft%ngm) - REAL(DP), INTENT(OUT) :: ga(3,dfft%nnr) - ! - INTEGER :: ipol - COMPLEX(DP), ALLOCATABLE :: gaux(:) - ! - ! - ALLOCATE( gaux( dfft%nnr ) ) - ! - ! ... multiply by (iG) to get (\grad_ipol a)(G) ... - ! - ga(:,:) = 0.D0 - ! - DO ipol = 1, 3 - ! - gaux(:) = (0.0_dp,0.0_dp) - ! - gaux(dfft%nl(:)) = g(ipol,:) * CMPLX( -AIMAG(a(:)), REAL(a(:)), kind=DP) - ! - IF ( gamma_only ) THEN - ! - gaux(dfft%nlm(:)) = CMPLX( REAL( gaux(dfft%nl(:)) ), & - -AIMAG( gaux(dfft%nl(:)) ), kind=DP) - ! - END IF - ! - ! ... bring back to R-space, (\grad_ipol a)(r) ... - ! - CALL invfft ('Rho', gaux, dfft) - ! - ! ...and add the factor 2\pi/a missing in the definition of G - ! - ga(ipol,:) = ga(ipol,:) + tpiba * REAL( gaux(:) ) - ! - END DO - ! - DEALLOCATE( gaux ) - ! - RETURN - ! -END SUBROUTINE gradrho -!---------------------------------------------------------------------------- -SUBROUTINE grad_dot( dfft, a, g, alat, da ) - !---------------------------------------------------------------------------- - ! - ! ... Calculates da = \sum_i \grad_i a_i in R-space - ! - USE cell_base, ONLY : tpiba - USE kinds, ONLY : DP - USE control_flags, ONLY : gamma_only - USE fft_interfaces,ONLY : fwfft, invfft - USE fft_types, ONLY : fft_type_descriptor - ! - IMPLICIT NONE - ! - TYPE(fft_type_descriptor),INTENT(IN) :: dfft - REAL(DP), INTENT(IN) :: a(3,dfft%nnr), g(3,dfft%ngm), alat - REAL(DP), INTENT(OUT) :: da(dfft%nnr) - ! - INTEGER :: n, ipol - COMPLEX(DP), ALLOCATABLE :: aux(:), gaux(:) - ! - ! - ALLOCATE( aux(dfft%nnr), gaux(dfft%nnr) ) - ! - gaux(:) = (0.0_dp,0.0_dp) - ! - DO ipol = 1, 3 - ! - aux = CMPLX( a(ipol,:), 0.0_dp, kind=DP) - ! - ! ... bring a(ipol,r) to G-space, a(G) ... - ! - CALL fwfft ('Rho', aux, dfft) - ! - DO n = 1, dfft%ngm - ! - gaux(dfft%nl(n)) = gaux(dfft%nl(n)) + g(ipol,n) * & - CMPLX( -AIMAG( aux(dfft%nl(n)) ), & - REAL( aux(dfft%nl(n)) ), kind=DP) - ! - END DO - ! - END DO - ! - IF ( gamma_only ) THEN - ! - DO n = 1, dfft%ngm - ! - gaux(dfft%nlm(n)) = CONJG( gaux(dfft%nl(n)) ) - ! - END DO - ! - END IF - ! - ! ... bring back to R-space, (\grad_ipol a)(r) ... - ! - CALL invfft ('Rho', gaux, dfft) - ! - ! ... add the factor 2\pi/a missing in the definition of G and sum - ! - da(:) = tpiba * REAL( gaux(:) ) - ! - DEALLOCATE( aux, gaux ) - ! - RETURN - ! -END SUBROUTINE grad_dot diff --git a/PW/src/loc_scdm.f90 b/PW/src/loc_scdm.f90 index 67731560c..c36e1918c 100644 --- a/PW/src/loc_scdm.f90 +++ b/PW/src/loc_scdm.f90 @@ -304,7 +304,7 @@ IMPLICIT NONE write(stdout,'(7x,A,f12.6)') 'DenMax = ', DenMax ! gradient on the exx grid - call fft_gradient( dfftt, den, gt, grad_den ) + call fft_gradient_r2r ( dfftt, den, gt, grad_den ) charge = Zero GrdAve = Zero GrdMax = Zero diff --git a/PW/src/stres_gradcorr.f90 b/PW/src/stres_gradcorr.f90 index f8273f9f4..90200fc5c 100644 --- a/PW/src/stres_gradcorr.f90 +++ b/PW/src/stres_gradcorr.f90 @@ -73,7 +73,7 @@ subroutine stres_gradcorr( rho, rhog, rho_core, rhog_core, kedtau, nspin, & rho(:,is) = fac * rho_core(:) + rho(:,is) rhog(:,is) = fac * rhog_core(:) + rhog(:,is) ! - CALL gradrho( dfft, rhog(1,is), g, grho(1,1,is) ) + CALL fft_gradient_g2r( dfft, rhog(1,is), g, grho(1,1,is) ) ! END DO ! diff --git a/PW/src/v_of_rho.f90 b/PW/src/v_of_rho.f90 index 4e25be531..cc6ffa068 100644 --- a/PW/src/v_of_rho.f90 +++ b/PW/src/v_of_rho.f90 @@ -111,7 +111,7 @@ SUBROUTINE v_xc_meta( rho, rho_core, rhog_core, etxc, vtxc, v, kedtaur ) USE fft_base, ONLY : dfftp USE gvect, ONLY : g, ngm USE lsda_mod, ONLY : nspin - USE cell_base, ONLY : omega, alat + USE cell_base, ONLY : omega USE spin_orb, ONLY : domag USE funct, ONLY : xc, xc_spin, tau_xc, tau_xc_spin, get_meta USE scf, ONLY : scf_type @@ -178,7 +178,7 @@ SUBROUTINE v_xc_meta( rho, rho_core, rhog_core, etxc, vtxc, v, kedtaur ) rhoout(:,is) = fac * rho_core(:) + rhoout(:,is) rhogsum(:,is) = fac * rhog_core(:) + rhogsum(:,is) ! - CALL gradrho( dfftp, rhogsum(1,is), g, grho(1,1,is) ) + CALL fft_gradient_g2r( dfftp, rhogsum(1,is), g, grho(1,1,is) ) ! END DO ! @@ -292,7 +292,7 @@ SUBROUTINE v_xc_meta( rho, rho_core, rhog_core, etxc, vtxc, v, kedtaur ) ! DO is = 1, nspin ! - CALL grad_dot( dfftp, h(1,1,is), g, alat, dh ) + CALL fft_graddot( dfftp, h(1,1,is), g, dh ) ! v(:,is) = v(:,is) - dh(:) !