diff --git a/GWW/simple/commutator.f90 b/GWW/simple/commutator.f90 index 33f52f466..dcada1a97 100644 --- a/GWW/simple/commutator.f90 +++ b/GWW/simple/commutator.f90 @@ -14,6 +14,7 @@ subroutine gen_beta_simple (qk, npw_max, dvkb) USE klist, ONLY : ngk USE gvect, ONLY : mill, eigts1, eigts2, eigts3, g USE uspp, ONLY : nkb, indv, nhtol, nhtolm + USE beta_mod, ONLY : interp_dbeta USE uspp_param, ONLY : lmaxkb, nbetam, nh USE io_global, ONLY : stdout ! @@ -137,6 +138,7 @@ subroutine gen_beta_simple_2 (qk, npw_max, u, dvkb) USE klist, ONLY : ngk, igk_k USE gvect, ONLY : mill, eigts1, eigts2, eigts3, g USE uspp, ONLY : nkb, indv, nhtol, nhtolm + USE beta_mod, ONLY : interp_beta USE uspp_param, ONLY : upf, lmaxkb, nbetam, nh ! implicit none diff --git a/PP/src/pw2gw.f90 b/PP/src/pw2gw.f90 index 3c48f6b6e..7a046bc14 100644 --- a/PP/src/pw2gw.f90 +++ b/PP/src/pw2gw.f90 @@ -150,7 +150,7 @@ SUBROUTINE compute_gw( omegamin, omegamax, d_omega, use_gmaps, qplda, vkb, vxcdi USE parallel_include USE scf, ONLY : rho, rho_core, rhog_core USE ener, ONLY : etxc, vtxc - + USE beta_mod, ONLY : interp_beta, interp_dbeta USE uspp_param, ONLY : upf, nh, nbetam USE ions_base, ONLY : ntyp => nsp USE klist, ONLY : ngk diff --git a/PW/src/hinit0.f90 b/PW/src/hinit0.f90 index 6b71f87e1..ca3b27c3c 100644 --- a/PW/src/hinit0.f90 +++ b/PW/src/hinit0.f90 @@ -22,6 +22,7 @@ SUBROUTINE hinit0() USE klist, ONLY : qnorm USE gvecw, ONLY : ecutwfc USE vlocal, ONLY : strf + USE beta_mod, ONLY : init_tab_beta USE realus, ONLY : generate_qpointlist, betapointlist, & init_realspace_vars, real_space USE ldaU, ONLY : lda_plus_U, Hubbard_projectors @@ -42,6 +43,7 @@ SUBROUTINE hinit0() ! IMPLICIT NONE REAL (dp) :: alat_old, qmax + INTEGER :: ierr LOGICAL :: is_tau_read = .FALSE. ! #if defined (__ENVIRON) @@ -60,10 +62,20 @@ SUBROUTINE hinit0() IF (tbeta_smoothing) CALL init_us_b0(ecutwfc,intra_bgrp_comm) IF (tq_smoothing) CALL init_us_0(ecutrho,intra_bgrp_comm) qmax = (qnorm + sqrt(ecutrho))*cell_factor - ! qmax is the maximum needed |q+G|, increased by a factor (20% or so) - ! to avoid too frequent reallocations in variable-cell calculations - ! (qnorm=max|q| may be needed for hybrid EXX or phonon calculations) + ! + ! qmax is the maximum |q+G|, for all G needed by the charge density, + ! increased by a factor (20% or so) to avoid too frequent reallocations + ! in variable-cell calculations ( norm is an estimate of max|q|, that + ! may be needed for hybrid EXX or phonon calculations) + ! CALL init_us_1(nat, ityp, omega, qmax, intra_bgrp_comm) + ! + ! fill interpolation table for beta functions + ! qmax as above, for all G needed by wavefunctions + ! + qmax = (qnorm + sqrt(ecutwfc))*cell_factor + CALL init_tab_beta ( qmax, omega, intra_bgrp_comm, ierr ) + ! IF ( lda_plus_U .AND. ( Hubbard_projectors == 'pseudo' ) ) CALL init_q_aeps() CALL init_tab_atwfc (omega, intra_bgrp_comm) ! diff --git a/PW/src/read_file_new.f90 b/PW/src/read_file_new.f90 index 60d6b2c9c..36ca5f6f8 100644 --- a/PW/src/read_file_new.f90 +++ b/PW/src/read_file_new.f90 @@ -282,11 +282,14 @@ SUBROUTINE post_xml_init ( ) USE rism_module, ONLY : rism_tobe_alive, rism_pot3d USE rism3d_facade, ONLY : lrism3d, rism3d_initialize, rism3d_read_to_restart USE xc_lib, ONLY : xclib_dft_is_libxc, xclib_init_libxc + USE beta_mod, ONLY : init_tab_beta + USE klist, ONLY : qnorm ! IMPLICIT NONE ! REAL(DP) :: ehart, etxc, vtxc, etotefield, charge, qmax CHARACTER(LEN=37) :: dft_name + INTEGER :: ierr ! ! ... initialize Libxc if needed ! @@ -359,9 +362,20 @@ SUBROUTINE post_xml_init ( ) CALL init_vloc() IF (tbeta_smoothing) CALL init_us_b0(ecutwfc,intra_bgrp_comm) IF (tq_smoothing) CALL init_us_0(ecutrho,intra_bgrp_comm) + ! + ! qmax is the maximum |G|, for all G needed by the charge density + ! qmax = sqrt(ecutrho)*cell_factor CALL init_us_1(nat, ityp, omega, qmax, intra_bgrp_comm) + ! + ! fill interpolation table for beta functions + ! qmax is the maximum |q+G|, for all G needed by the wavefunctions + ! + qmax = (qnorm + sqrt(ecutwfc))*cell_factor + CALL init_tab_beta ( qmax, omega, intra_bgrp_comm, ierr ) + ! IF ( lda_plus_u .AND. ( Hubbard_projectors == 'pseudo' ) ) CALL init_q_aeps() + ! CALL init_tab_atwfc(omega, intra_bgrp_comm) ! CALL struc_fact( nat, tau, nsp, ityp, ngm, g, bg, dfftp%nr1, dfftp%nr2,& diff --git a/PW/src/scale_h.f90 b/PW/src/scale_h.f90 index fa8ae6431..ea599871d 100644 --- a/PW/src/scale_h.f90 +++ b/PW/src/scale_h.f90 @@ -29,6 +29,7 @@ SUBROUTINE scale_h USE rism_module, ONLY : lrism, rism_reinit3d USE mp, ONLY : mp_max USE mp_bands, ONLY : intra_bgrp_comm + USE beta_mod, ONLY : scale_tab_beta USE vloc_mod, ONLY : scale_tab_vloc USE rhoc_mod, ONLY : scale_tab_rhc USE rhoat_mod, ONLY : scale_tab_rhoat @@ -86,6 +87,7 @@ SUBROUTINE scale_h ! scale the non-local pseudopotential tables ! call scale_uspp_data( omega_old/omega ) + call scale_tab_beta( omega_old/omega ) CALL scale_tab_rhc( omega_old/omega ) CALL scale_tab_rhoat( omega_old/omega ) CALL scale_tab_qrad( omega_old/omega ) diff --git a/upflib/CMakeLists.txt b/upflib/CMakeLists.txt index 2e68c5b64..e12a2f524 100644 --- a/upflib/CMakeLists.txt +++ b/upflib/CMakeLists.txt @@ -8,7 +8,7 @@ set(src_upflib init_us_b0.f90 init_us_1.f90 init_tab_atwfc.f90 - init_tab_beta.f90 + beta_mod.f90 interp_atwfc.f90 init_us_2_acc.f90 gth.f90 diff --git a/upflib/Makefile b/upflib/Makefile index 1abb76748..98bba3e2a 100644 --- a/upflib/Makefile +++ b/upflib/Makefile @@ -10,10 +10,11 @@ MODFLAGS= $(MOD_FLAG)../UtilXlib $(MOD_FLAG). # OBJS_GPU are GPU-specific OBJS_DEP= \ -vloc_mod.o \ -rhoc_mod.o \ -rhoat_mod.o \ +beta_mod.o \ qrad_mod.o \ +rhoat_mod.o \ +rhoc_mod.o \ +vloc_mod.o \ qvan2.o \ dqvan2.o \ gen_us_dj.o \ @@ -22,8 +23,7 @@ init_us_0.o \ init_us_b0.o \ init_us_1.o \ init_us_2_acc.o \ -init_tab_atwfc.o \ -init_tab_beta.o +init_tab_atwfc.o OBJS_NODEP= \ atom.o \ diff --git a/upflib/beta_mod.f90 b/upflib/beta_mod.f90 new file mode 100644 index 000000000..13d5f545a --- /dev/null +++ b/upflib/beta_mod.f90 @@ -0,0 +1,236 @@ +! +! Copyright (C) 2024 Quantum ESPRESSO Foundation +! This file is distributed under the terms of the +! GNU General Public License. See the file `License' +! in the root directory of the present distribution, +! or http://www.gnu.org/copyleft/gpl.txt . +! +! +MODULE beta_mod + ! + !! Variables and routines for nonlocal beta functions in numerical form + !! Contains generation of interpolation tables in reciprocal space, + !! interpolation routines and other utility routines + !! Code moved to upflib and restructured by Paolo Giannozzi, 2024 + ! + USE upf_kinds, ONLY : dp + USE upf_const, ONLY : fpi, e2 + ! + IMPLICIT NONE + ! + PRIVATE + PUBLIC :: init_tab_beta + PUBLIC :: deallocate_tab_beta + PUBLIC :: scale_tab_beta + PUBLIC :: interp_beta + PUBLIC :: interp_dbeta + ! + SAVE + ! + INTEGER :: nqx = 0 + !! size of interpolation table + REAL(DP), PARAMETER:: dq = 0.01_dp + !! grid step for interpolation table + REAL(DP) :: qmax = 0.0_dp + !! max q covered by the interpolation table + REAL(DP), ALLOCATABLE :: tab_beta(:,:,:) + !! interpolation table for numerical beta functions in reciprocal space + ! +CONTAINS +! +!---------------------------------------------------------------------- +SUBROUTINE init_tab_beta ( qmax_, omega, comm, ierr ) + !---------------------------------------------------------------------- + ! + ! Compute interpolation table for beta(G) radial functions + ! + USE upf_kinds, ONLY : dp + USE atom, ONLY : rgrid + USE uspp_param, ONLY : upf, lmaxq, nbetam, nsp + USE mp, ONLY : mp_sum + USE m_gth, ONLY : mk_ffnl_gth + ! + IMPLICIT NONE + ! + REAL(dp), INTENT(IN) :: qmax_ + !! Interpolate q up to qmax_ (sqrt(Ry), q^2 is an energy) + REAL(dp), INTENT(IN) :: omega + !! Unit-cell volume + INTEGER, INTENT(IN) :: comm + !! MPI communicator, to split the workload + INTEGER, INTENT(OUT) :: ierr + !! error code: ierr = 0 if interpolation table (IT) was allocated + !! ierr =-1 if IT had insufficient size and was re-allocated + !! ierr =-2 if IT was already present and nothing is done + ! + INTEGER :: ndm, startq, lastq, nt, l, nb, iq, ir + REAL(dp) :: qi + ! q-point grid for interpolation + REAL(dp) :: pref + ! the prefactor of the Q functions + real(DP) :: vqint, d1 + ! + REAL(dp), allocatable :: aux (:) + ! work space + REAL(dp), allocatable :: besr(:) + ! work space + ! + IF ( .NOT. ALLOCATED(tab_beta) ) THEN + !! table not yet allocated + qmax = qmax_ + ierr = 0 + ELSE IF ( qmax_ > qmax ) THEN + !! table ìs allocated but dimension insufficient: re-allocate + !! (with some margin so that this does not happen too often) + !$acc exit data delete(tab_beta) + DEALLOCATE ( tab_beta ) + qmax = qmax_ + MAX(dq*100,qmax_-qmax) + ierr =-1 + ELSE + !! table already computed: exit + ierr =-2 + RETURN + END IF + nqx = INT( qmax/dq + 4) + allocate(tab_beta(nqx,nbetam,nsp)) + !$acc enter data create(tab_beta) + ndm = MAXVAL ( upf(:)%kkbeta ) + allocate( aux (ndm) ) + allocate (besr( ndm)) + pref = fpi / sqrt (omega) + call divide (comm, nqx, startq, lastq) + tab_beta (:,:,:) = 0.d0 + do nt = 1, nsp + do nb = 1, upf(nt)%nbeta + l = upf(nt)%lll (nb) + do iq = startq, lastq + qi = (iq - 1) * dq + if ( upf(nt)%is_gth ) then + CALL mk_ffnl_gth( nt, nb, 1, omega, [ qi ] , tab_beta(iq,nb,nt) ) + else + call sph_bes (upf(nt)%kkbeta, rgrid(nt)%r, qi, l, besr) + do ir = 1, upf(nt)%kkbeta + aux (ir) = upf(nt)%beta (ir, nb) * besr (ir) * rgrid(nt)%r(ir) + enddo + call simpson (upf(nt)%kkbeta, aux, rgrid(nt)%rab, vqint) + tab_beta (iq, nb, nt) = vqint * pref + end if + enddo + enddo + enddo + deallocate (besr) + deallocate (aux) + ! + call mp_sum( tab_beta, comm ) +!$acc update device (tab_beta) + ! +END SUBROUTINE init_tab_beta +! +!---------------------------------------------------------------------- +SUBROUTINE interp_beta( nt, npw_, qg, vq ) + !---------------------------------------------------------------------- + ! + USE upf_kinds, ONLY : dp + USE uspp_param, ONLY : upf, nbetam + ! + implicit none + integer, intent(in) :: nt, npw_ + real(dp), intent(in ) :: qg(npw_) + real(dp), intent(out) :: vq(npw_,nbetam) + ! + integer :: i0, i1, i2, i3, nbnt, nb, ig + real(dp):: qgr, px, ux, vx, wx + ! + nbnt = upf(nt)%nbeta + !$acc data present (tab_beta) present_or_copyin (qg) present_or_copyout (vq) + !$acc parallel loop collapse(2) + do nb = 1, nbnt + DO ig = 1, npw_ + qgr = qg(ig) + px = qgr / dq - DBLE(INT(qgr/dq)) + ux = 1.0_dp - px + vx = 2.0_dp - px + wx = 3.0_dp - px + i0 = INT(qgr/dq) + 1 + i1 = i0 + 1 + i2 = i0 + 2 + i3 = i0 + 3 + if ( i3 <= nqx ) then + vq(ig,nb) = & + tab_beta(i0,nb,nt) * ux * vx * wx / 6.0_dp + & + tab_beta(i1,nb,nt) * px * vx * wx / 2.0_dp - & + tab_beta(i2,nb,nt) * px * ux * wx / 2.0_dp + & + tab_beta(i3,nb,nt) * px * ux * vx / 6.0_dp + else + !! This case should never happen if tab_beta is properly allocated + !! (setting q_max to be large enough) - for compatibility with GWW + vq(ig,nb) = 0.0_dp + end if + END DO + END DO + !$acc end data + !---------------------------------------------------------------------- +END SUBROUTINE interp_beta +!---------------------------------------------------------------------- +! +!---------------------------------------------------------------------- +SUBROUTINE interp_dbeta( nt, npw, qg, vq ) + !---------------------------------------------------------------------- + ! + USE upf_kinds, ONLY : dp + USE uspp_param, ONLY : upf, nbetam + ! + implicit none + integer, intent(in) :: nt, npw + real(dp), intent(in ) :: qg(npw) + real(dp), intent(out) :: vq(npw,nbetam) + ! + integer :: i0, i1, i2, i3, nbnt, nb, ig + real(dp):: qgr, px, ux, vx, wx + ! + nbnt = upf(nt)%nbeta + !$acc data present (tab_beta) present_or_copyin (qg) present_or_copyout (vq) + !$acc parallel loop collapse(2) + DO nb = 1, nbnt + DO ig = 1, npw + qgr = qg(ig) + px = qgr / dq - INT(qgr/dq) + ux = 1.0_dp - px + vx = 2.0_dp - px + wx = 3.0_dp - px + i0 = qgr / dq + 1 + i1 = i0 + 1 + i2 = i0 + 2 + i3 = i0 + 3 + IF ( i3 <= nqx ) THEN + vq(ig,nb) = ( tab_beta(i0,nb,nt) * (-vx*wx-ux*wx-ux*vx)/6.0_dp + & + tab_beta(i1,nb,nt) * (+vx*wx-px*wx-px*vx)/2.0_dp - & + tab_beta(i2,nb,nt) * (+ux*wx-px*wx-px*ux)/2.0_dp + & + tab_beta(i3,nb,nt) * (+ux*vx-px*vx-px*ux)/6.0_dp )/dq + ELSE + vq(ig,nb) = 0.0_dp + END IF + ENDDO + END DO + !$acc end data + !---------------------------------------------------------------------- +END SUBROUTINE interp_dbeta +! + subroutine deallocate_tab_beta () + implicit none + ! + !$acc exit data delete(tab_beta) + if( allocated( tab_beta ) ) deallocate( tab_beta ) + ! + end subroutine deallocate_tab_beta + ! + subroutine scale_tab_beta( vol_ratio_m1 ) + ! vol_ratio_m1 = omega_old / omega + implicit none + real(DP), intent(in) :: vol_ratio_m1 + ! + tab_beta(:,:,:) = tab_beta(:,:,:) * SQRT(vol_ratio_m1) +!$acc update device ( tab_beta) + end subroutine scale_tab_beta + ! +END MODULE beta_mod diff --git a/upflib/gen_us_dj.f90 b/upflib/gen_us_dj.f90 index ef09b51c1..b73892547 100644 --- a/upflib/gen_us_dj.f90 +++ b/upflib/gen_us_dj.f90 @@ -18,6 +18,7 @@ SUBROUTINE gen_us_dj_base( npw, npwx, igk, xk, nat, tau, ityp, ntyp, tpiba, & USE upf_const, ONLY: tpi USE uspp, ONLY: nkb, indv, nhtol, nhtolm USE uspp_param, ONLY: lmaxkb, nbetam, nh, nhm + USE beta_mod, ONLY: interp_dbeta ! IMPLICIT NONE ! @@ -202,47 +203,3 @@ SUBROUTINE gen_us_dj_base( npw, npwx, igk, xk, nat, tau, ityp, ntyp, tpiba, & RETURN ! END SUBROUTINE gen_us_dj_base -! -!---------------------------------------------------------------------- -SUBROUTINE interp_dbeta( nt, npw, qg, vq ) - !---------------------------------------------------------------------- - ! - USE upf_kinds, ONLY : dp - USE uspp_param, ONLY : upf, nbetam - USE uspp_data, ONLY : nqx, dq, tab_beta - ! - implicit none - integer, intent(in) :: nt, npw - real(dp), intent(in ) :: qg(npw) - real(dp), intent(out) :: vq(npw,nbetam) - ! - integer :: i0, i1, i2, i3, nbnt, nb, ig - real(dp):: qgr, px, ux, vx, wx - ! - nbnt = upf(nt)%nbeta - !$acc data present (tab_beta) present_or_copyin (qg) present_or_copyout (vq) - !$acc parallel loop collapse(2) - DO nb = 1, nbnt - DO ig = 1, npw - qgr = qg(ig) - px = qgr / dq - INT(qgr/dq) - ux = 1.0_dp - px - vx = 2.0_dp - px - wx = 3.0_dp - px - i0 = qgr / dq + 1 - i1 = i0 + 1 - i2 = i0 + 2 - i3 = i0 + 3 - IF ( i3 <= nqx ) THEN - vq(ig,nb) = ( tab_beta(i0,nb,nt) * (-vx*wx-ux*wx-ux*vx)/6.0_dp + & - tab_beta(i1,nb,nt) * (+vx*wx-px*wx-px*vx)/2.0_dp - & - tab_beta(i2,nb,nt) * (+ux*wx-px*wx-px*ux)/2.0_dp + & - tab_beta(i3,nb,nt) * (+ux*vx-px*vx-px*ux)/6.0_dp )/dq - ELSE - vq(ig,nb) = 0.0_dp - END IF - ENDDO - END DO - !$acc end data - !---------------------------------------------------------------------- -END SUBROUTINE interp_dbeta diff --git a/upflib/gen_us_dy.f90 b/upflib/gen_us_dy.f90 index 4dea2bfc6..c7f8f666e 100644 --- a/upflib/gen_us_dy.f90 +++ b/upflib/gen_us_dy.f90 @@ -20,6 +20,7 @@ SUBROUTINE gen_us_dy_base( npw, npwx, igk, xk, nat, tau, ityp, ntyp, tpiba, & USE upf_const, ONLY: tpi USE uspp, ONLY: nkb, indv, nhtol, nhtolm USE uspp_param, ONLY: upf, lmaxkb, nbetam, nh, nhm + USE beta_mod, ONLY: interp_beta ! IMPLICIT NONE ! diff --git a/upflib/init_tab_beta.f90 b/upflib/init_tab_beta.f90 deleted file mode 100644 index b34253297..000000000 --- a/upflib/init_tab_beta.f90 +++ /dev/null @@ -1,69 +0,0 @@ -! -! Copyright (C) 2021-2023 Quantum ESPRESSO Foundation -! This file is distributed under the terms of the -! GNU General Public License. See the file `License' -! in the root directory of the present distribution, -! or http://www.gnu.org/copyleft/gpl.txt . -! -!---------------------------------------------------------------------- -SUBROUTINE init_tab_beta ( omega, intra_bgrp_comm ) - !---------------------------------------------------------------------- - ! - ! Compute interpolation table for beta(G) radial functions - ! - USE upf_kinds, ONLY : dp - USE upf_const, ONLY : fpi - USE atom, ONLY : rgrid - USE uspp_param, ONLY : upf, lmaxq, nbetam, nsp - USE uspp_data, ONLY : nqx, dq, tab_beta - USE mp, ONLY : mp_sum - USE m_gth, ONLY : mk_ffnl_gth - ! - IMPLICIT NONE - ! - real(DP), intent(in) :: omega - integer, intent(in) :: intra_bgrp_comm - ! - INTEGER :: ndm, startq, lastq, nt, l, nb, iq, ir - REAL(dp) :: qi - ! q-point grid for interpolation - REAL(dp) :: pref - ! the prefactor of the Q functions - real(DP) :: vqint, d1 - ! - REAL(dp), allocatable :: aux (:) - ! work space - REAL(dp), allocatable :: besr(:) - ! work space - ! - ndm = MAXVAL ( upf(:)%kkbeta ) - allocate( aux (ndm) ) - allocate (besr( ndm)) - pref = fpi / sqrt (omega) - call divide (intra_bgrp_comm, nqx, startq, lastq) - tab_beta (:,:,:) = 0.d0 - do nt = 1, nsp - do nb = 1, upf(nt)%nbeta - l = upf(nt)%lll (nb) - do iq = startq, lastq - qi = (iq - 1) * dq - if ( upf(nt)%is_gth ) then - CALL mk_ffnl_gth( nt, nb, 1, omega, [ qi ] , tab_beta(iq,nb,nt) ) - else - call sph_bes (upf(nt)%kkbeta, rgrid(nt)%r, qi, l, besr) - do ir = 1, upf(nt)%kkbeta - aux (ir) = upf(nt)%beta (ir, nb) * besr (ir) * rgrid(nt)%r(ir) - enddo - call simpson (upf(nt)%kkbeta, aux, rgrid(nt)%rab, vqint) - tab_beta (iq, nb, nt) = vqint * pref - end if - enddo - enddo - enddo - deallocate (besr) - deallocate (aux) - ! - call mp_sum( tab_beta, intra_bgrp_comm ) -!$acc update device (tab_beta) - ! -END SUBROUTINE init_tab_beta diff --git a/upflib/init_us_1.f90 b/upflib/init_us_1.f90 index 0fa4947dc..db3d8959c 100644 --- a/upflib/init_us_1.f90 +++ b/upflib/init_us_1.f90 @@ -251,10 +251,6 @@ subroutine init_us_1( nat, ityp, omega, qmax, intra_bgrp_comm ) end do end if ! - ! fill interpolation table for beta functions - ! - CALL init_tab_beta ( omega, intra_bgrp_comm ) - ! #if defined __CUDA ! ! update GPU memory (taking care of zero-dim allocations) diff --git a/upflib/init_us_2_acc.f90 b/upflib/init_us_2_acc.f90 index 3f5157041..1cec50c87 100644 --- a/upflib/init_us_2_acc.f90 +++ b/upflib/init_us_2_acc.f90 @@ -20,6 +20,7 @@ SUBROUTINE init_us_2_acc( npw_, npwx, igk_, q_, nat, tau, ityp, & USE upf_const, ONLY : tpi USE uspp, ONLY : nkb, nhtol, nhtolm, indv USE uspp_param, ONLY : lmaxkb, nbetam, nhm, nh, nsp + USE beta_mod, ONLY : interp_beta ! implicit none ! @@ -168,52 +169,6 @@ SUBROUTINE init_us_2_acc( npw_, npwx, igk_, q_, nat, tau, ityp, & deallocate(vkb1) ! return + !---------------------------------------------------------------------- end subroutine init_us_2_acc -! -!---------------------------------------------------------------------- -SUBROUTINE interp_beta( nt, npw_, qg, vq ) - !---------------------------------------------------------------------- - ! - USE upf_kinds, ONLY : dp - USE uspp_param, ONLY : upf, nbetam - USE uspp_data, ONLY : nqx, dq, tab_beta - ! - implicit none - integer, intent(in) :: nt, npw_ - real(dp), intent(in ) :: qg(npw_) - real(dp), intent(out) :: vq(npw_,nbetam) - ! - integer :: i0, i1, i2, i3, nbnt, nb, ig - real(dp):: qgr, px, ux, vx, wx - ! - nbnt = upf(nt)%nbeta - !$acc data present (tab_beta) present_or_copyin (qg) present_or_copyout (vq) - !$acc parallel loop collapse(2) - do nb = 1, nbnt - DO ig = 1, npw_ - qgr = qg(ig) - px = qgr / dq - DBLE(INT(qgr/dq)) - ux = 1.0_dp - px - vx = 2.0_dp - px - wx = 3.0_dp - px - i0 = INT(qgr/dq) + 1 - i1 = i0 + 1 - i2 = i0 + 2 - i3 = i0 + 3 - if ( i3 <= nqx ) then - vq(ig,nb) = & - tab_beta(i0,nb,nt) * ux * vx * wx / 6.0_dp + & - tab_beta(i1,nb,nt) * px * vx * wx / 2.0_dp - & - tab_beta(i2,nb,nt) * px * ux * wx / 2.0_dp + & - tab_beta(i3,nb,nt) * px * ux * vx / 6.0_dp - else - !! This case should never happen if tab_beta is properly allocated - !! (setting q_max to be large enough) - for compatibility with GWW - vq(ig,nb) = 0.0_dp - end if - END DO - END DO - !$acc end data - !---------------------------------------------------------------------- -END SUBROUTINE interp_beta !---------------------------------------------------------------------- diff --git a/upflib/uspp_data.f90 b/upflib/uspp_data.f90 index a01701dc7..d2daf7269 100644 --- a/upflib/uspp_data.f90 +++ b/upflib/uspp_data.f90 @@ -16,7 +16,7 @@ MODULE uspp_data PRIVATE ! PUBLIC :: nqxq, nqx, dq - PUBLIC :: tab_beta, tab_at + PUBLIC :: tab_at ! PUBLIC :: allocate_uspp_data PUBLIC :: deallocate_uspp_data @@ -27,9 +27,7 @@ MODULE uspp_data INTEGER :: nqx !! number of interpolation points REAL(DP), PARAMETER:: dq = 0.01D0 - !! space between points in the pseudopotential tab_beta. - REAL(DP), ALLOCATABLE :: tab_beta(:,:,:) - !! interpolation table for PP projectorss + !! space between interpolation points REAL(DP), ALLOCATABLE :: tab_at(:,:,:) !! interpolation table for atomic wfc ! @@ -43,8 +41,6 @@ contains if (nqxq_/=nqxq) call upf_error("allocate_uspp_data","invalid nqxq_",1) if (nqx_/=nqx) call upf_error("allocate_uspp_data","invalid nqx_",1) ! - allocate(tab_beta(nqx_,nbetam,nsp)) - !$acc enter data create(tab_beta) allocate(tab_at(nqx_,nwfcm,nsp)) !$acc enter data create(tab_at) ! @@ -52,8 +48,6 @@ contains ! subroutine deallocate_uspp_data() implicit none - !$acc exit data delete(tab_beta) - if( allocated( tab_beta ) ) deallocate( tab_beta ) !$acc exit data delete(tab_at) if( allocated( tab_at ) ) deallocate( tab_at ) ! @@ -64,9 +58,8 @@ contains implicit none real(DP), intent(in) :: vol_ratio_m1 ! - tab_beta(:,:,:) = tab_beta(:,:,:) * SQRT(vol_ratio_m1) tab_at(:,:,:) = tab_at(:,:,:) * SQRT(vol_ratio_m1) -!$acc update device (tab_at, tab_beta) +!$acc update device (tab_at) end subroutine scale_uspp_data ! END MODULE uspp_data