Commit 818d55e4 authored by Samuel Poncé's avatar Samuel Poncé

USPP support for EPW

Ultrasoft support for EPW + various cleaning.
Addition of a new uspp test.

- Roxana Margine and Samuel Ponce
parent 4c453cb2
......@@ -37,6 +37,7 @@ transport.o \
transport_iter.o \
wigner.o \
a2f.o \
adddvscf2.o \
allocate_epwq.o \
bcast_epw_input.o \
bloch2wan.o \
......@@ -67,6 +68,7 @@ loadqmesh.o \
loadumat.o \
load_rebal.o \
nesting_fn.o \
newdq2.o \
openfilepw.o \
rgd_blk_epw_fine_mem.o \
plot_band.o \
......
!
! Copyright (C) 2010-2016 Samuel Ponce', Roxana Margine, Carla Verdi, Feliciano Giustino
!
! Copyright (C) 2001-2003 PWSCF group
! 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 .
!
! adapted from LR_Modules/adddvscf.f90 (QE)
!
!----------------------------------------------------------------------
SUBROUTINE adddvscf2( ipert, ik )
!----------------------------------------------------------------------
!!
!! This routine computes the contribution of the selfconsistent
!! change of the potential to the known part of the linear
!! system and adds it to dvpsi.
!! It implements the second term in Eq. B30 of PRB 64, 235118 (2001).
!! Only used in the case of USPP.
!!
!! Roxana Margine - Jan 2019: Updated based on QE 6.3
!! SP - Jan 2019: Clean
!!
USE kinds, ONLY : DP
USE uspp_param, ONLY : upf, nh
USE uspp, ONLY : vkb, okvan
USE lsda_mod, ONLY : lsda, current_spin, isk
USE ions_base, ONLY : ntyp => nsp, nat, ityp
USE wvfct, ONLY : npwx
USE lrus, ONLY : int3, int3_nc, becp1
USE qpoint, ONLY : npwq
USE eqv, ONLY : dvpsi
USE elph2, ONLY : lower_band, upper_band
USE noncollin_module, ONLY : noncolin, npol
USE constants_epw, ONLY : czero
!
IMPLICIT NONE
!
INTEGER, INTENT(in) :: ik
!! Counter on k-point
INTEGER, INTENT(in) :: ipert
!! Counter on Vscf perturbations
!
! Local variables
!
INTEGER :: na
!! Counter on atoms
INTEGER :: nt
!! Counter on atomic types
INTEGER :: ibnd
!! Counter on bands
INTEGER :: ih
!! Counter on beta functions
INTEGER :: jh
!! Counter on beta functions
INTEGER :: ijkb0
!! Auxiliary variable for counting
INTEGER :: ikb
!! Counter on becp functions
INTEGER :: jkb
!! Counter on becp functions
INTEGER :: is
!! Counter on polarization
INTEGER :: js
!! Counter on polarization
INTEGER :: ijs
!! Counter on combined is and js polarization
!
COMPLEX(DP) :: sumA
!! auxiliary sum variable
COMPLEX(DP) :: sum_nc(npol)
!! auxiliary sum variable non-collinear case
!
IF (.not.okvan) RETURN
!
CALL start_clock('adddvscf2')
!
IF (lsda) current_spin = isk(ik)
!
ijkb0 = 0
DO nt = 1, ntyp
IF ( upf(nt)%tvanp ) THEN
DO na = 1, nat
IF (ityp(na) .eq. nt) THEN
!
! we multiply the integral for the becp term and the beta_n
!
DO ibnd = lower_band, upper_band
do ih = 1, nh(nt)
ikb = ijkb0 + ih
IF (noncolin) THEN
sum_nc = czero
ELSE
sumA = czero
ENDIF
DO jh = 1, nh(nt)
jkb = ijkb0 + jh
IF (noncolin) THEN
ijs = 0
DO is = 1, npol
DO js = 1, npol
ijs = ijs + 1
sum_nc(is) = sum_nc(is) + int3_nc(ih,jh,na,ijs,ipert) * &
becp1(ik)%nc(jkb,js,ibnd)
ENDDO
ENDDO
ELSE
sumA = sumA + int3(ih,jh,na,current_spin,ipert) * &
becp1(ik)%k(jkb,ibnd)
ENDIF
ENDDO
IF (noncolin) THEN
CALL zaxpy( npwq, sum_nc(1), vkb(1,ikb), 1, dvpsi(1,ibnd), 1 )
CALL zaxpy( npwq, sum_nc(2), vkb(1,ikb), 1, dvpsi(1+npwx,ibnd), 1 )
ELSE
CALL zaxpy( npwq, sumA, vkb(1,ikb), 1, dvpsi(1,ibnd), 1 )
ENDIF
ENDDO
ENDDO
ijkb0 = ijkb0 + nh(nt)
ENDIF
ENDDO
ELSE
DO na = 1, nat
IF (ityp(na) .eq. nt) ijkb0 = ijkb0 + nh(nt)
ENDDO
ENDIF
ENDDO
!
CALL stop_clock('adddvscf2')
!
RETURN
!
END SUBROUTINE adddvscf2
......@@ -17,8 +17,8 @@
!! response problem
!!
!! RM - Nov/Dec - 2014 - Imported the noncolinear case implemented by xlzhang
!! SP - 2016 - Update for QE 5
!! RM - Dec 2018 - Updated based on QE 6.3
!! SP - 2016 - Updated for QE 5
!! RM - Jan 2019 - Updated based on QE 6.3
!!
USE ions_base, ONLY : nat, ntyp => nsp
USE pwcom, ONLY : npwx, nbnd, nspin
......@@ -36,7 +36,6 @@
USE becmod, ONLY : becp, allocate_bec_type
USE uspp_param, ONLY : nhm
USE uspp, ONLY : okvan, nkb
USE units_ph, ONLY : this_dvkb3_is_on_file, this_pcxpsi_is_on_file
USE modes, ONLY : u, npert, name_rap_mode, num_rap_mode
USE klist, ONLY : nks
USE fft_base, ONLY : dfftp
......@@ -45,8 +44,6 @@
!
IMPLICIT NONE
!
! SP: Had to add these allocations becaue they are now private in QE 5.0.
! See private 'becp_nc' variable from espresso/Modules/becmod.f90.
INTEGER :: ik
!! k-point
INTEGER :: ipol
......@@ -56,7 +53,7 @@
!
ALLOCATE (evq(npwx*npol, nbnd))
ALLOCATE (dpsi ( npwx*npol, nbnd))
ALLOCATE( transp_temp(nstemp) )
ALLOCATE (transp_temp(nstemp))
!
ALLOCATE (vlocq(ngm, ntyp))
! SP: nrxx is not used in QE 5 ==> tg_nnr is the maximum among nnr
......@@ -69,7 +66,8 @@
! nnr_tg = local number of grid elements for task group FFT ( ~nr1*nr2*nr3/proc3 )
! --> tg = task group
! ALLOCATE (dmuxc ( dffts%nnr, nspin, nspin))
ALLOCATE (dmuxc( dfftp%nnr, nspin_mag, nspin_mag))
!
ALLOCATE (dmuxc(dfftp%nnr, nspin_mag, nspin_mag))
ALLOCATE (eigqts(nat))
ALLOCATE (rtau(3, 48, nat))
ALLOCATE (u(3 * nat, 3 * nat))
......@@ -82,31 +80,27 @@
ALLOCATE (int4(nhm * (nhm + 1)/2, 3, 3, nat, nspin_mag))
ALLOCATE (int5(nhm * (nhm + 1)/2, 3, 3, nat , nat))
IF (noncolin) THEN
ALLOCATE(int1_nc(nhm, nhm, 3, nat, nspin))
ALLOCATE(int4_nc(nhm, nhm, 3, 3, nat, nspin))
ALLOCATE(becsum_nc(nhm*(nhm+1)/2, nat, npol, npol))
ALLOCATE(alphasum_nc(nhm*(nhm+1)/2, 3, nat, npol, npol))
ALLOCATE (int1_nc(nhm, nhm, 3, nat, nspin))
ALLOCATE (int4_nc(nhm, nhm, 3, 3, nat, nspin))
ALLOCATE (becsum_nc(nhm*(nhm+1)/2, nat, npol, npol))
ALLOCATE (alphasum_nc(nhm*(nhm+1)/2, 3, nat, npol, npol))
IF (lspinorb) THEN
ALLOCATE(int2_so(nhm, nhm, 3, nat, nat, nspin))
ALLOCATE(int5_so(nhm, nhm, 3, 3, nat, nat, nspin))
ALLOCATE (int2_so(nhm, nhm, 3, nat, nat, nspin))
ALLOCATE (int5_so(nhm, nhm, 3, 3, nat, nat, nspin))
ENDIF
ENDIF ! noncolin
ALLOCATE (alphasum (nhm * (nhm + 1)/2, 3, nat , nspin_mag))
ALLOCATE (this_dvkb3_is_on_file(nks))
this_dvkb3_is_on_file(:) = .false.
ALLOCATE (alphasum(nhm * (nhm + 1)/2, 3, nat, nspin_mag))
ENDIF
ALLOCATE (this_pcxpsi_is_on_file(nks,3))
this_pcxpsi_is_on_file(:,:) = .false.
!
ALLOCATE (becp1(nks))
ALLOCATE (alphap(3,nks))
!
DO ik = 1, nks
CALL allocate_bec_type(nkb, nbnd, becp1(ik))
DO ipol = 1, 3
CALL allocate_bec_type(nkb, nbnd, alphap(ipol,ik))
ENDDO
END DO
CALL allocate_bec_type(nkb, nbnd, becp1(ik))
DO ipol = 1, 3
CALL allocate_bec_type(nkb, nbnd, alphap(ipol,ik))
ENDDO
ENDDO
CALL allocate_bec_type(nkb, nbnd, becp)
!
IF (elph) ALLOCATE (el_ph_mat(nbnd, nbnd, nks, 3*nat))
......
......@@ -91,9 +91,9 @@
!! $$\mathbf{r}\cdot\mathbf{k}
REAL(kind=DP) :: tmp
!! Maximum value of the real space Hamiltonian
REAL(kind=dp) :: et_opt(nbnd,nks)
REAL(kind=DP) :: et_opt(nbnd,nks)
!! hamiltonian eigenvalues within the outer window in the first ndimwin(ik) entries
REAL(kind=dp) :: et_tmp(nbnd,nks)
REAL(kind=DP) :: et_tmp(nbnd,nks)
!! temporary array for hamiltonian eigenvalues
!
COMPLEX(kind=DP) :: chs(nbndsub, nbndsub, nks)
......@@ -264,7 +264,7 @@
USE mp_world, ONLY : mpime
USE mp, ONLY : mp_barrier, mp_sum
!
implicit none
IMPLICIT NONE
!
! input variables
!
......@@ -459,7 +459,7 @@
! Check spatial decay of Dipole in Wannier basis
! the unit in r-space is angstrom
!
IF (mpime.eq.ionode_id) then
IF (mpime.eq.ionode_id) THEN
OPEN(unit=iudecayP,file='decay.P')
WRITE(iudecayP, '(/3x,a/)') '#Spatial decay of dipole in Wannier basis'
DO ir = 1, nrr
......@@ -515,7 +515,7 @@
USE mp, ONLY : mp_barrier
USE mp_global, ONLY : inter_pool_comm
!
implicit none
IMPLICIT NONE
!
! input variables
!
......@@ -636,7 +636,9 @@
USE mp_global, ONLY : inter_pool_comm, my_pool_id
USE mp, ONLY : mp_barrier, mp_sum
USE mp_world, ONLY : mpime
implicit none
USE division, ONLY : fkbounds
!
IMPLICIT NONE
!
! input variables
!
......@@ -750,7 +752,7 @@
!
! setup rotation matrix - we need access to all for the k+b
cu_big = czero
CALL ckbounds(ikstart, ikstop)
CALL fkbounds(nkstot, ikstart, ikstop)
cu_big(:,:,ikstart:ikstop) = cu(:,:,:)
CALL mp_sum(cu_big,inter_pool_comm)
!
......@@ -1050,7 +1052,8 @@
USE mp_global, ONLY : inter_pool_comm
USE mp , ONLY : mp_sum
USE mp_world, ONLY : mpime
implicit none
!
IMPLICIT NONE
!
! input variables
!
......@@ -1213,7 +1216,7 @@
USE mp, ONLY : mp_barrier
USE mp_world, ONLY : mpime
!
implicit none
IMPLICIT NONE
!
! Input variables
!
......@@ -1338,7 +1341,8 @@
USE io_global, ONLY : ionode_id
USE mp, ONLY : mp_barrier
USE mp_world, ONLY : mpime
implicit none
!
IMPLICIT NONE
!
! input variables - note irvec is dimensioned with nrr_k
! (which is assumed to be larger than nrr_q)
......@@ -1359,7 +1363,7 @@
!! Coordinates of real space vector
!
REAL(kind=DP), INTENT(in) :: xk(3, nq)
!! Kpoint coordinates (cartesian in units of 2piba)
!! K-point coordinates (cartesian in units of 2piba)
!
COMPLEX(kind=DP), INTENT(in) :: epmatwe (nbnd, nbnd, nrr_k, nmodes)
!! EP matrix in electron-wannier representation and phonon bloch representation
......
......@@ -13,16 +13,7 @@
!! This subroutine is called from elphon_shuffle_wrap for each
!! nq1*nq2*nq3 phonon on the coarse mesh.
!!
!! obsolete:
!! variables for folding of the k+q mesh into the kmesh
!! there may be at most 9 different G_0 to refold k+q into k
!! of the first BZ (including G_0=0)
!! [actually there are 3^3=27 translations but for a fixed q only
!! 2^3=8 out of these 27 are possible since q has a definite direction]
!!
!! new:
!! I need all q-vectors in the same run, so I consider all the
!! 125 possibilities by setting ng0vec = 125 in ../PW/set_kplusq.f90
!! It folds the k+q mesh into the k mesh using 5^3 G_0 translations
!!
!! SP - 2016 - iverbosity cannot be tested here. Generates Tb of data ...
!
......@@ -49,7 +40,7 @@
INTEGER :: ik
!! k-point index
INTEGER :: jk
!! other k-point index
!! k-point index
INTEGER :: i
!! Index of the k+q on the k-grid
INTEGER :: j
......@@ -103,10 +94,6 @@
abs(zz-nint(zz)) .le. eps5
IF (.not.in_the_list) CALL errore('createkmap','q-vec not commensurate',1)
!
! previously the 27 possible translations
! now we have 125 possibilities for LSCO, because the unit cell is far from cubic... sob!
! ng0vec must be redefined
!
ng0vec = 0
DO ig1 = -2, 2
DO ig2 = -2, 2
......@@ -310,9 +297,9 @@
xx_c(ik) = xk(1,ik) * nk1
yy_c(ik) = xk(2,ik) * nk2
zz_c(ik) = xk(3,ik) * nk3
in_the_list = abs(xx-nint(xx)).le.eps5 .and. &
abs(yy-nint(yy)).le.eps5 .and. &
abs(zz-nint(zz)).le.eps5
in_the_list = abs(xx_c(ik)-nint(xx_c(ik))) .le. eps5 .AND. &
abs(yy_c(ik)-nint(yy_c(ik))) .le. eps5 .AND. &
abs(zz_c(ik)-nint(zz_c(ik))) .le. eps5
IF (.not.in_the_list) CALL errore('createkmap2','is this a uniform k-mesh?',1)
!
IF ( (xx_c(ik) .lt. -eps5) .OR. (yy_c(ik) .lt. -eps5) .OR. (zz_c(ik) .lt. -eps5) ) &
......@@ -328,9 +315,9 @@
xx = xkq(1,ik) * nk1
yy = xkq(2,ik) * nk2
zz = xkq(3,ik) * nk3
in_the_list = abs(xx-nint(xx)).le.eps5 .and. &
abs(yy-nint(yy)).le.eps5 .and. &
abs(zz-nint(zz)).le.eps5
in_the_list = abs(xx-nint(xx)) .le. eps5 .AND. &
abs(yy-nint(yy)) .le. eps5 .AND. &
abs(zz-nint(zz)) .le. eps5
IF (.not.in_the_list) CALL errore('createkmap2','k+q does not fall on k-grid',1)
!
! find the index of this k+q in the k-grid
......@@ -433,16 +420,16 @@
!! itoj(i) returns the index of the G-vector in the sorted list
!! that was at i-th position in the unsorted list
!
REAL(KIND=DP) :: xx, yy, zz
!
REAL(DP), ALLOCATABLE :: gg(:)
REAL(kind=DP) :: xx, yy, zz
!! k-point in crystal coords. in multiple of nk1, nk2, nk3
REAL(kind=DP), ALLOCATABLE :: gg(:)
!! G^2 in increasing order
REAL(DP), ALLOCATABLE :: g(:,:)
REAL(kind=DP), ALLOCATABLE :: g(:,:)
!! G-vectors cartesian components in increasing order of G^2
REAL(DP), ALLOCATABLE :: g2sort_g(:)
REAL(kind=DP), ALLOCATABLE :: g2sort_g(:)
!! G-vectors for the current processor
REAL(DP) :: tx(3), ty(3), t(3)
REAL(DP), ALLOCATABLE :: tt(:)
REAL(kind=DP) :: tx(3), ty(3), t(3)
REAL(kind=DP), ALLOCATABLE :: tt(:)
!
LOGICAL :: in_the_list
!! Variable in the list
......@@ -638,7 +625,7 @@
CALL mp_barrier(inter_pool_comm)
CALL mp_barrier(inter_image_comm)
!
DEALLOCATE( ig_l2g, mill, igsrt, jtoi, itoj, g, gg )
DEALLOCATE( ig_l2g, mill, mill_unsorted, igsrt, jtoi, itoj, g, gg )
!
RETURN
!
......@@ -693,11 +680,10 @@
INTEGER :: ig1_use, ig2_use, ig2_guess, notfound, guess_skip
!! Temporary G-vectors indices
INTEGER :: indold, indnew
!! Counter on G_0 vectors Indices for writing to file
!! Counter on G_0 vectors indices for writing to file
!
LOGICAL :: tfound
!
!
ALLOCATE( gmap(ngm_g,ng0vec) )
gmap(:,:) = 0
guess_skip = 0
......@@ -714,8 +700,6 @@
IF (indnew.ne.indold) WRITE(stdout,'(a)',advance='no') '#'
indold = indnew
!
!
!
notfound = 0
DO ig1 = 1, ngm_g
!
......@@ -786,11 +770,6 @@
IF (iukgmap .ne. stdout) CLOSE(iukgmap)
WRITE(stdout,*)
!
end subroutine refold
RETURN
!
END SUBROUTINE refold
......@@ -17,7 +17,7 @@
!! This routine calculates two integrals of the Q functions and
!! its derivatives with c V_loc and V_eff which are used
!! to compute term dV_bare/dtau * psi in addusdvqpsi.
!! The result is stored in int1,int2. The routine is called
!! The result is stored in int1, int2, int4, int5. The routine is called
!! for each q in nqc.
!!
!! RM - Nov/Dec 2014
......@@ -188,7 +188,7 @@
IF (jpol >= ipol) THEN
DO ig = 1, ngm
aux3(ig) = aux5(ig) * &
( g(jpol,ig) + xq(jpol) )
( g(jpol,ig) + xq(jpol) )
ENDDO
int5(ijh,ipol,jpol,na,nb) = &
CONJG(fact) * tpiba2 * omega * &
......@@ -253,8 +253,8 @@
ENDDO
CALL mp_sum(int1, intra_pool_comm)
CALL mp_sum(int2, intra_pool_comm)
call mp_sum(int4, intra_pool_comm)
call mp_sum(int5, intra_pool_comm)
CALL mp_sum(int4, intra_pool_comm)
CALL mp_sum(int5, intra_pool_comm)
IF (noncolin) THEN
CALL set_int12_nc(0)
int4_nc = czero
......@@ -292,4 +292,4 @@
CALL stop_clock ('dvanqq2')
RETURN
!
END SUBROUTINE dvanqq2
END SUBROUTINE dvanqq2
This diff is collapsed.
This diff is collapsed.
This diff is collapsed.
......@@ -7,60 +7,95 @@
! present distribution, or http://www.gnu.org/copyleft.gpl.txt .
!
!-----------------------------------------------------------------------
SUBROUTINE elphon_shuffle ( iq_irr, nqc_irr, iq, gmapsym, eigv, isym, xq0, timerev )
SUBROUTINE elphon_shuffle( iq_irr, nqc_irr, iq, gmapsym, eigv, isym, xq0, timerev )
!-----------------------------------------------------------------------
!!
!! Electron-phonon calculation from data saved in fildvscf
!! Shuffle2 mode (shuffle on electrons + load all phonon q's)
!!
!! no ultrasoft yet
!! No ultrasoft yet
!!
!! RM - Nov/Dec 2014
!! Imported the noncolinear case implemented by xlzhang
!!
!
!! Roxana Margine - Jan 2019: Updated based on QE 6.3 for US
!!
!-----------------------------------------------------------------------
!
USE kinds, ONLY : DP
USE mp, ONLY : mp_barrier, mp_sum
USE mp_global, ONLY : my_pool_id, nproc_pool,npool,kunit,&
inter_pool_comm
USE mp_global, ONLY : my_pool_id, nproc_pool, npool, kunit, &
inter_pool_comm
USE mp_images, ONLY : nproc_image
USE ions_base, ONLY : nat
USE pwcom, ONLY : nbnd, nks, nkstot
USE gvect, ONLY : ngm
USE gvecs, ONLY : doublegrid
USE kinds, ONLY : DP
USE modes, ONLY : nmodes, nirr, npert, u
USE elph2, ONLY : epmatq, el_ph_mat
USE constants_epw, ONLY : czero, cone
USE lrus, ONLY : int3, int3_nc
USE uspp, ONLY : okvan
USE lsda_mod, ONLY : nspin
USE fft_base, ONLY : dfftp, dffts
USE fft_interfaces, ONLY : fft_interpolate
USE noncollin_module, ONLY : nspin_mag
! USE noncollin_module, ONLY : noncolin
!
implicit none
!
integer :: irr, imode0, ipert, is, iq, iq_irr, nqc_irr
! counter on the representations
! counter on the modes
! the change of Vscf due to perturbations
! the current qpoint in the uniform grid
! the current ireducible qpoint
! the total number of irreducible qpoints in the list
complex(kind=DP), POINTER :: dvscfin(:,:,:), dvscfins (:,:,:)
logical :: timerev
! true if we are using time reversal
!
integer :: tmp_pool_id, ik0, ik, ibnd, jbnd
integer :: iks, nkl, nkr
integer :: gmapsym ( ngm, 48 ), isym
! the correspondence G-->S(G)
! the symmetry which generates the current q in the star
complex(kind=DP) :: eigv (ngm, 48)
! e^{ iGv} for 1...nsym ( v the fractional translation)
real(kind=DP) :: xq0(3)
! the first q in the star (cartesian)
!
CALL start_clock ('elphon')
USE uspp_param, ONLY : nhm
USE constants_epw, ONLY : czero, cone
USE fft_interfaces, ONLY : fft_interpolate
USE noncollin_module, ONLY : nspin_mag, noncolin
!
IMPLICIT NONE
!
INTEGER, INTENT(in) :: iq_irr
!! Current ireducible q-point
INTEGER, INTENT(in) :: nqc_irr
!! Total number of irreducible q-points in the list
INTEGER, INTENT(in) :: iq
!! Current q-point in the star of iq_irr q-point
INTEGER, INTENT(in) :: gmapsym(ngm,48)
!! Correspondence G-->S(G)
INTEGER, INTENT(in) :: isym
!! The symmetry which generates the current q in the star
REAL(DP), INTENT(in) :: xq0(3)
!! The first q-point in the star (cartesian coords.)
COMPLEX(DP), INTENT(in) :: eigv(ngm,48)
!! e^{iGv} for 1...nsym (v the fractional translation)
LOGICAL, INTENT(in) :: timerev
!! true if we are using time reversal
!
! Local variables
!
INTEGER :: irr
!! Counter on representations
INTEGER :: imode0
!! Counter on modes
INTEGER :: ipert
!! Change of Vscf due to perturbations
INTEGER :: npe
!! Number of perturbations for irr representation
INTEGER :: is
!! Counter on spin
INTEGER :: ik
!! Counter on k-points in the pool
INTEGER :: ibnd
!! Counter on bands
INTEGER :: jbnd
!! Counter on bands
INTEGER :: tmp_pool_id
!! temporary pool id
INTEGER :: iks
!! Index of the first k-point block in this pool
INTEGER :: ik0
!! Index of iks - 1
INTEGER :: nkl
!!
INTEGER :: nkr
!!
!
COMPLEX(DP), POINTER :: dvscfin(:,:,:)
!! Change of the scf potential
COMPLEX(DP), POINTER :: dvscfins(:,:,:)
!! Change of the scf potential (smooth)
!
CALL start_clock('elphon_shuffle')
!
ik0 = 0
tmp_pool_id = 0
......@@ -69,7 +104,7 @@
IF (npool.gt.1) THEN
!
! number of kpoint blocks, kpoints per pool and reminder
kunit = 1
kunit = 1
nkl = kunit * ( nkstot / npool )
nkr = ( nkstot - nkl * npool ) / kunit
! the reminder goes to the first nkr pools
......@@ -89,40 +124,44 @@
!
imode0 = 0
DO irr = 1, nirr
ALLOCATE (dvscfin ( dfftp%nnr , nspin_mag , npert(irr)) )
!DBSP
! if (noncolin) then
! ALLOCATE (dvscfin ( dfftp%nnr , nspin_mag , npert(irr)) )
! endif
!END
npe = npert(irr)
ALLOCATE( dvscfin(dfftp%nnr, nspin_mag, npe) )
IF (okvan) THEN
ALLOCATE( int3(nhm, nhm, nat, nspin_mag, npe) )
IF (noncolin) ALLOCATE( int3_nc(nhm, nhm, nat, nspin, npe) )
ENDIF
!
! read the <prefix>.dvscf_q[iq] files
!
dvscfin = (0.d0,0.d0)
IF (my_pool_id.eq.0) THEN
DO ipert = 1, npert (irr)
CALL readdvscf ( dvscfin(1,1,ipert), imode0 + ipert, iq_irr, nqc_irr )
dvscfin = czero
IF ( my_pool_id.eq.0 ) THEN
DO ipert = 1, npe
CALL readdvscf( dvscfin(1,1,ipert), imode0 + ipert, iq_irr, nqc_irr )
ENDDO
ENDIF
CALL mp_sum(dvscfin,inter_pool_comm)
!
!
IF (doublegrid) THEN
ALLOCATE (dvscfins ( dffts%nnr , nspin_mag , npert(irr)) )
DO is = 1, nspin_mag
DO ipert = 1, npert(irr)
CALL fft_interpolate (dfftp, dvscfin(:,is,ipert), dffts, dvscfins(:,is,ipert))
ENDDO
ENDDO