mirror of https://gitlab.com/QEF/q-e.git
Final cft3/cft3s => fwfft/invfft conversion
git-svn-id: http://qeforge.qe-forge.org/svn/q-e/trunk/espresso@7048 c92efa57-630b-4861-b058-cf58834340f0
This commit is contained in:
parent
31bb050a09
commit
f2f3cae392
|
@ -60,7 +60,7 @@ all : tldeps tddfpt.x
|
|||
#
|
||||
tddfpt.x : $(LROBJS) $(PWOBJS) $(QEMODS) $(PHOBJS) $(PWOBJS) $(EEOBJS) $(LIBOBJS)
|
||||
$(LD) $(LDFLAGS) -o $@ \
|
||||
$(LROBJS) $(PWOBJS) $(QEMODS) $(PHOBJS) $(EEOBJS) $(LIBOBJS) $(LIBS)
|
||||
$(LROBJS) $(PHOBJS) $(PWOBJS) $(EEOBJS) $(QEMODS) $(LIBOBJS) $(LIBS)
|
||||
- ( cd ../../bin ; ln -fs ../TDDFPT/src/$@ . )
|
||||
- if [ -d ../bin ] ; then ( cd ../bin ; ln -fs ../src/$@ . ); fi
|
||||
|
||||
|
|
|
@ -5,12 +5,12 @@ subroutine lr_apply_liouvillian( evc1, evc1_new, sevc1_new, interaction )
|
|||
! OBM: or to be more exact this function is responsible for calculating L.q(i) and (L^T).p(i)
|
||||
! where q is evc1 and return partial qdash(i+1) and pdash(i+1) in evc1_new . Ultrasoft additions are
|
||||
! handled here...
|
||||
! interaction=.true. corresponds to eq.(32)
|
||||
! interaction=.true. corresponds to eq.(32)
|
||||
! interaction=.false. corresponds to eq.(33)
|
||||
! in Ralph Gebauer, Brent Walker J. Chem. Phys., 127, 164106 (2007)
|
||||
!---------------------------------------------------------------------
|
||||
! in Ralph Gebauer, Brent Walker J. Chem. Phys., 127, 164106 (2007)
|
||||
!---------------------------------------------------------------------
|
||||
!
|
||||
! Modified by Osman Baris Malcioglu in 2009
|
||||
! Modified by Osman Baris Malcioglu in 2009
|
||||
#include "f_defs.h"
|
||||
!
|
||||
use ions_base, only : ityp, nat, ntyp=>nsp
|
||||
|
@ -28,7 +28,7 @@ subroutine lr_apply_liouvillian( evc1, evc1_new, sevc1_new, interaction )
|
|||
use uspp, only : vkb, nkb, okvan
|
||||
use uspp_param, only : nhm, nh
|
||||
use wavefunctions_module, only : psic
|
||||
use wvfct, only : nbnd, npwx, igk, g2kin, et
|
||||
use wvfct, only : nbnd, npwx, igk, g2kin, et
|
||||
use control_flags, only : gamma_only
|
||||
USE realus, ONLY : real_space, fft_orbital_gamma, initialisation_level, &
|
||||
bfft_orbital_gamma, calbec_rs_gamma, add_vuspsir_gamma, &
|
||||
|
@ -38,7 +38,7 @@ subroutine lr_apply_liouvillian( evc1, evc1_new, sevc1_new, interaction )
|
|||
USE io_global, ONLY : stdout
|
||||
USE dfunct, ONLY : newq
|
||||
USE control_flags, ONLY : tqr
|
||||
|
||||
|
||||
|
||||
!
|
||||
implicit none
|
||||
|
@ -88,16 +88,16 @@ subroutine lr_apply_liouvillian( evc1, evc1_new, sevc1_new, interaction )
|
|||
dvrs(:,1)=rho_1(:,1)
|
||||
!
|
||||
!call lr_dv_of_drho(dvrs)
|
||||
allocate( dvrs_temp(nrxx, nspin) )
|
||||
allocate( dvrs_temp(nrxx, nspin) )
|
||||
dvrs_temp=CMPLX(dvrs,0.0d0) !OBM: This memory copy was hidden in lr_dv_of_drho, can it be avoided?
|
||||
call dv_of_drho(0,dvrs_temp,.false.)
|
||||
dvrs=DBLE(dvrs_temp)
|
||||
deallocate(dvrs_temp)
|
||||
!
|
||||
if ( okvan ) then
|
||||
if ( okvan ) then
|
||||
if ( tqr ) then
|
||||
call newq_r(dvrs,d_deeq,.true.)
|
||||
else
|
||||
else
|
||||
call newq(dvrs,d_deeq,.true.)
|
||||
endif
|
||||
endif
|
||||
|
@ -135,7 +135,7 @@ subroutine lr_apply_liouvillian( evc1, evc1_new, sevc1_new, interaction )
|
|||
!
|
||||
! Tamm-dancoff interaction
|
||||
!
|
||||
write(stdout,'(5X,"lr_apply_liouvillian: applying interaction: tamm-dancoff")')
|
||||
write(stdout,'(5X,"lr_apply_liouvillian: applying interaction: tamm-dancoff")')
|
||||
!
|
||||
! Here evc1_new contains the interaction
|
||||
!
|
||||
|
@ -217,7 +217,7 @@ contains
|
|||
if ( interaction ) then
|
||||
!
|
||||
call start_clock('interaction')
|
||||
|
||||
|
||||
if (nkb > 0 .and. okvan) then
|
||||
! calculation of becp2
|
||||
becp2(:,:) = 0.0d0
|
||||
|
@ -245,7 +245,7 @@ contains
|
|||
enddo
|
||||
!
|
||||
enddo
|
||||
!
|
||||
!
|
||||
enddo
|
||||
!
|
||||
ijkb0 = ijkb0 + nh(nt)
|
||||
|
@ -300,8 +300,8 @@ contains
|
|||
DO jh = 1, nh(nt)
|
||||
!
|
||||
jkb = ijkb0 + jh
|
||||
w1(ih) = w1(ih) + becp2(jkb, ibnd)
|
||||
IF ( ibnd+1 .le. nbnd ) w2(ih) = w2(ih) + becp2(jkb, ibnd+1)
|
||||
w1(ih) = w1(ih) + becp2(jkb, ibnd)
|
||||
IF ( ibnd+1 .le. nbnd ) w2(ih) = w2(ih) + becp2(jkb, ibnd+1)
|
||||
!
|
||||
END DO
|
||||
!
|
||||
|
@ -357,9 +357,9 @@ contains
|
|||
!
|
||||
call h_psi(npwx,npw_k(1),nbnd,evc1(1,1,1),sevc1_new(1,1,1))
|
||||
!
|
||||
! spsi1 = s*evc1
|
||||
! spsi1 = s*evc1
|
||||
!
|
||||
if (real_space_debug > 9 ) then
|
||||
if (real_space_debug > 9 ) then
|
||||
do ibnd=1,nbnd,2
|
||||
call fft_orbital_gamma(evc1(:,:,1),ibnd,nbnd)
|
||||
call s_psir_gamma(ibnd,nbnd)
|
||||
|
@ -405,7 +405,7 @@ contains
|
|||
call start_clock('interaction')
|
||||
!
|
||||
! evc1_new is used as a container for the interaction
|
||||
!
|
||||
!
|
||||
evc1_new(:,:,:)=(0.0d0,0.0d0)
|
||||
!
|
||||
do ik=1,nks
|
||||
|
@ -422,7 +422,7 @@ contains
|
|||
!
|
||||
! Back to reciprocal space
|
||||
!
|
||||
call cft3s(psic,nr1s,nr2s,nr3s,nrx1s,nrx2s,nrx3s,-2)
|
||||
CALL fwfft ('Wave', psic, dffts)
|
||||
!
|
||||
do ig=1,npw_k(ik)
|
||||
!
|
||||
|
@ -478,7 +478,7 @@ contains
|
|||
!
|
||||
enddo
|
||||
!
|
||||
!evc1_new(ik) = evc1_new(ik) + vkb*becp2(ik)
|
||||
!evc1_new(ik) = evc1_new(ik) + vkb*becp2(ik)
|
||||
call zgemm( 'N', 'N', npw_k(ik), nbnd, nkb, (1.d0,0.d0), vkb, &
|
||||
npwx, becp2, nkb, (1.d0,0.d0), evc1_new(:,:,ik), npwx )
|
||||
!
|
||||
|
|
|
@ -1,13 +1,13 @@
|
|||
!-----------------------------------------------------------------------
|
||||
subroutine lr_calc_dens( evc1, response_calc )
|
||||
!---------------------------------------------------------------------
|
||||
! ... calculates response charge density from linear response
|
||||
! ... calculates response charge density from linear response
|
||||
! ... orbitals and ground state orbitals
|
||||
!---------------------------------------------------------------------
|
||||
!
|
||||
! Modified by Osman Baris Malcioglu in 2009
|
||||
! Modified by Osman Baris Malcioglu in 2009
|
||||
!
|
||||
! Input : evc1 (qdash etc) Output: rho_1 (=2*sum_v (revc0_v(r) . revc1_v(r,w) v:valance state index, r denotes a transformation to real space)
|
||||
! Input : evc1 (qdash etc) Output: rho_1 (=2*sum_v (revc0_v(r) . revc1_v(r,w) v:valance state index, r denotes a transformation t
|
||||
! in case of uspps becsum is also calculated here
|
||||
! in case of charge response calculation, the rho_tot is calculated here
|
||||
|
||||
|
@ -54,14 +54,14 @@ subroutine lr_calc_dens( evc1, response_calc )
|
|||
logical, intent(in) :: response_calc
|
||||
!
|
||||
! functions
|
||||
real(kind=dp) :: ddot
|
||||
real(kind=dp) :: ddot
|
||||
!
|
||||
! local variables
|
||||
integer :: ir,ik,ibnd,jbnd,ig,ijkb0,np,na,ijh,ih,jh,ikb,jkb,ispin
|
||||
integer :: i, j, k, l
|
||||
integer :: i, j, k, l
|
||||
real(kind=dp) :: w1,w2,scal
|
||||
real(kind=dp) :: rho_sum!,weight
|
||||
real(kind=dp), allocatable :: rho_sum_resp_x(:),rho_sum_resp_y(:),rho_sum_resp_z(:) ! These are temporary buffers for response charge storage
|
||||
real(kind=dp), allocatable :: rho_sum_resp_x(:),rho_sum_resp_y(:),rho_sum_resp_z(:) ! These are temporary buffers for response cha
|
||||
|
||||
!complex(kind=dp), external :: ZDOTC
|
||||
!complex(kind=dp), allocatable :: spsi(:,:)
|
||||
|
@ -88,9 +88,9 @@ subroutine lr_calc_dens( evc1, response_calc )
|
|||
else
|
||||
call lr_calc_dens_k()
|
||||
endif
|
||||
!print *, "rho_1 after lr_calc_dens calculates",SUM(rho_1)
|
||||
!print *, "rho_1 after lr_calc_dens calculates",SUM(rho_1)
|
||||
!print *, "norm of evc1 after lr_calc_dens calculates", lr_dot(evc1(1,1,1),evc1(1,1,1))
|
||||
|
||||
|
||||
!
|
||||
! ... If a double grid is used, interpolate onto the fine grid
|
||||
!
|
||||
|
@ -99,18 +99,18 @@ subroutine lr_calc_dens( evc1, response_calc )
|
|||
! ... Here we add the Ultrasoft contribution to the charge
|
||||
!
|
||||
!IF ( okvan ) CALL lr_addusdens(rho_1)
|
||||
!print *, "rho_1 before addusdens",SUM(rho_1)
|
||||
!print *, "rho_1 before addusdens",SUM(rho_1)
|
||||
!call start_clock('lrcd_usdens') !TQR makes a huge gain here
|
||||
if(okvan) then
|
||||
if (tqr) then
|
||||
CALL addusdens_r(rho_1,.false.)
|
||||
else
|
||||
CALL addusdens(rho_1)
|
||||
CALL addusdens(rho_1)
|
||||
endif
|
||||
endif
|
||||
!call stop_clock('lrcd_usdens')
|
||||
!
|
||||
!print *, "rho_1 after addusdens",SUM(rho_1)
|
||||
!print *, "rho_1 after addusdens",SUM(rho_1)
|
||||
#ifdef __PARA
|
||||
!call poolreduce(nrxx,rho_1)
|
||||
call mp_sum(rho_1, inter_pool_comm)
|
||||
|
@ -119,7 +119,7 @@ subroutine lr_calc_dens( evc1, response_calc )
|
|||
! check response charge density sums to 0
|
||||
!call start_clock('lrcd_sp') !Minimal lag, no need to improve
|
||||
if (lr_verbosity > 0) then
|
||||
|
||||
|
||||
do ispin = 1, nspin_mag
|
||||
rho_sum=0.0d0
|
||||
rho_sum=SUM(rho_1(:,ispin))
|
||||
|
@ -182,7 +182,7 @@ endif
|
|||
#endif
|
||||
write(stdout,'(5X,"Dumping plane sums of densities for iteration ",I4)') LR_iteration
|
||||
!
|
||||
filename = trim(prefix) // ".density_x"
|
||||
filename = trim(prefix) // ".density_x"
|
||||
tempfile = trim(tmp_dir) // trim(filename)
|
||||
!
|
||||
open (158, file = tempfile, form = 'formatted', status = 'unknown', position = 'append')
|
||||
|
@ -227,9 +227,9 @@ endif
|
|||
IF (charge_response == 2 .and. response_calc) then
|
||||
if (LR_iteration < itermax) WRITE(stdout,'(5x,"Calculating total response charge density")')
|
||||
! the charge response, it is actually equivalent to an element of
|
||||
! V^T . phi_v where V^T is the is the transpose of the Krylov subspace generated
|
||||
! V^T . phi_v where V^T is the is the transpose of the Krylov subspace generated
|
||||
! by the Lanczos algorithm. The total charge density can be written
|
||||
! as
|
||||
! as
|
||||
! \sum_(lanczos iterations) (V^T.phi_v) . w_T
|
||||
! Where w_T is the corresponding eigenvector from the solution of
|
||||
! (w-L)e_1 = w_T
|
||||
|
@ -254,7 +254,7 @@ endif
|
|||
call daxpy(nrxx, dble(w_T(LR_iteration)),rho_1(:,1),1,rho_1_tot(:,1),1) !spin not implemented
|
||||
endif
|
||||
If (lr_verbosity > 9) THEN
|
||||
if (LR_iteration == 2) then
|
||||
if (LR_iteration == 2) then
|
||||
call lr_dump_rho_tot_cube(rho_1(:,1),"first-rho1")
|
||||
endif
|
||||
if (LR_iteration == itermax .or. LR_iteration == itermax-1) call lr_dump_rho_tot_cube(rho_1(:,1),"last--rho1")
|
||||
|
@ -285,7 +285,7 @@ contains
|
|||
real_space_debug
|
||||
|
||||
!
|
||||
|
||||
|
||||
do ibnd=1,nbnd,2
|
||||
call fft_orbital_gamma(evc1(:,:,1),ibnd,nbnd)
|
||||
!
|
||||
|
@ -298,30 +298,30 @@ contains
|
|||
endif
|
||||
!call start_clock('lrcd-lp1')
|
||||
! OBM:
|
||||
! (n'(r,w)=2*sum_v (psi_v(r) . q_v(r,w))
|
||||
! (n'(r,w)=2*sum_v (psi_v(r) . q_v(r,w))
|
||||
! where psi are the ground state valance orbitals
|
||||
! and q_v are the standart batch representation (rotated)
|
||||
! and q_v are the standart batch representation (rotated)
|
||||
! response orbitals
|
||||
! Here, since the ith iteration is the best approximate we have
|
||||
! for the most dominant eigenvalues/vectors, an estimate for the response
|
||||
! charge density can be calculated. This is in no way the final
|
||||
! response charge density.
|
||||
! the loop is over real space points.
|
||||
! response charge density.
|
||||
! the loop is over real space points.
|
||||
do ir=1,nrxxs
|
||||
rho_1(ir,1)=rho_1(ir,1) &
|
||||
+2.0d0*(w1*real(revc0(ir,ibnd,1),dp)*real(psic(ir),dp)&
|
||||
+w2*aimag(revc0(ir,ibnd,1))*aimag(psic(ir)))
|
||||
enddo
|
||||
!
|
||||
!
|
||||
!call stop_clock('lrcd-lp1')
|
||||
! OBM - psic now contains the response functions at
|
||||
! OBM - psic now contains the response functions at
|
||||
! real space, eagerly putting all the real space stuff at this point.
|
||||
! notice that betapointlist() is called in lr_readin at the very start
|
||||
IF ( real_space_debug > 6 .and. okvan) then
|
||||
! The rbecp term
|
||||
call calbec_rs_gamma(ibnd,nbnd,becp%r)
|
||||
endif
|
||||
! End of real space stuff
|
||||
endif
|
||||
! End of real space stuff
|
||||
enddo
|
||||
!
|
||||
! ... If we have a US pseudopotential we compute here the becsum term
|
||||
|
@ -359,7 +359,7 @@ contains
|
|||
do ih = 1, nh(np)
|
||||
!
|
||||
ikb = ijkb0 + ih
|
||||
!
|
||||
!
|
||||
becsum(ijh,na,current_spin) = &
|
||||
becsum(ijh,na,current_spin) + &
|
||||
2.d0 * w1 * becp%r(ikb,ibnd) * becp1(ikb,ibnd)
|
||||
|
@ -428,7 +428,7 @@ contains
|
|||
psic(nls(igk_k(ig,ik)))=evc1(ig,ibnd,ik)
|
||||
enddo
|
||||
!
|
||||
call cft3s(psic,nr1s,nr2s,nr3s,nrx1s,nrx2s,nrx3s,2)
|
||||
CALL invfft ('Wave', psic, dffts)
|
||||
!
|
||||
w1=wg(ibnd,ik)/omega
|
||||
!
|
||||
|
@ -453,7 +453,7 @@ contains
|
|||
becsum(:,:,:) = 0.0d0
|
||||
!
|
||||
IF ( nkb > 0 .and. okvan ) then
|
||||
! call ccalbec(nkb,npwx,npw_k(ik),nbnd,becp,vkb,evc1)
|
||||
! call ccalbec(nkb,npwx,npw_k(ik),nbnd,becp,vkb,evc1)
|
||||
call calbec(npw_k(ik),vkb,evc1(:,:,ik),becp)
|
||||
endif
|
||||
!
|
||||
|
|
|
@ -55,7 +55,7 @@ subroutine lr_h_psiq (lda, n, m, psi, hpsi, spsi)
|
|||
!OBM debug
|
||||
!real(DP) :: obm_debug
|
||||
!complex(kind=dp), external :: ZDOTC
|
||||
|
||||
|
||||
|
||||
call start_clock ('h_psiq')
|
||||
If (lr_verbosity > 5) WRITE(stdout,'("<lr_h_psiq>")')
|
||||
|
@ -66,7 +66,7 @@ subroutine lr_h_psiq (lda, n, m, psi, hpsi, spsi)
|
|||
endif
|
||||
call stop_clock ('h_psiq')
|
||||
return
|
||||
contains
|
||||
contains
|
||||
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
|
||||
!k point part
|
||||
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
|
||||
|
@ -107,14 +107,14 @@ contains
|
|||
psic_nc(nls(igk(j)),1) = psi (j, ibnd)
|
||||
psic_nc(nls(igk(j)),2) = psi (j+lda, ibnd)
|
||||
enddo
|
||||
call cft3s (psic_nc(1,1), nr1s, nr2s, nr3s, nrx1s, nrx2s, nrx3s, 2)
|
||||
call cft3s (psic_nc(1,2), nr1s, nr2s, nr3s, nrx1s, nrx2s, nrx3s, 2)
|
||||
CALL invfft ('Wave', psic_nc(1,1), dffts)
|
||||
CALL invfft ('Wave', psic_nc(1,2), dffts)
|
||||
ELSE
|
||||
psic(:) = (0.d0, 0.d0)
|
||||
do j = 1, n
|
||||
psic (nls(igk(j))) = psi (j, ibnd)
|
||||
enddo
|
||||
call cft3s (psic, nr1s, nr2s, nr3s, nrx1s, nrx2s, nrx3s, 2)
|
||||
CALL invfft ('Wave', psic, dffts)
|
||||
END IF
|
||||
call stop_clock ('firstfft')
|
||||
!
|
||||
|
@ -146,8 +146,8 @@ contains
|
|||
!
|
||||
call start_clock ('secondfft')
|
||||
IF (noncolin) THEN
|
||||
call cft3s(psic_nc(1,1),nr1s,nr2s,nr3s,nrx1s,nrx2s,nrx3s,-2)
|
||||
call cft3s(psic_nc(1,2),nr1s,nr2s,nr3s,nrx1s,nrx2s,nrx3s,-2)
|
||||
CALL fwfft ('Wave', psic_nc(1,1), dffts)
|
||||
CALL fwfft ('Wave', psic_nc(1,2), dffts)
|
||||
!
|
||||
! addition to the total product
|
||||
!
|
||||
|
@ -156,7 +156,7 @@ contains
|
|||
hpsi (j+lda, ibnd) = hpsi (j+lda, ibnd) + psic_nc (nls(igk(j)), 2)
|
||||
enddo
|
||||
ELSE
|
||||
call cft3s (psic, nr1s, nr2s, nr3s, nrx1s, nrx2s, nrx3s, - 2)
|
||||
CALL fwfft ('Wave', psic, dffts)
|
||||
!
|
||||
! addition to the total product
|
||||
!
|
||||
|
@ -188,13 +188,13 @@ contains
|
|||
use uspp, only : nkb
|
||||
|
||||
IMPLICIT NONE
|
||||
|
||||
|
||||
call start_clock ('init')
|
||||
!
|
||||
! Here we apply the kinetic energy (k+G)^2 psi
|
||||
!
|
||||
if(gstart==2) psi(1,:)=cmplx(real(psi(1,:),dp),0.0d0,dp)
|
||||
!
|
||||
!
|
||||
!!OBM debug
|
||||
! obm_debug=0
|
||||
! do ibnd=1,m
|
||||
|
@ -209,7 +209,7 @@ contains
|
|||
do j=1,n
|
||||
hpsi(j,ibnd)=g2kin(j)*psi(j,ibnd)
|
||||
enddo
|
||||
enddo
|
||||
enddo
|
||||
!!OBM debug
|
||||
! obm_debug=0
|
||||
! do ibnd=1,m
|
||||
|
@ -260,9 +260,9 @@ contains
|
|||
! print *, "lr_h_psiq becp", obm_debug
|
||||
!!obm_debug
|
||||
|
||||
|
||||
|
||||
call add_vuspsi (lda, n, m, hpsi)
|
||||
!END IF
|
||||
!END IF
|
||||
!!OBM debug
|
||||
! obm_debug=0
|
||||
! do ibnd=1,m
|
||||
|
|
|
@ -1,7 +1,7 @@
|
|||
!-----------------------------------------------------------------------
|
||||
subroutine lr_read_wf()
|
||||
!---------------------------------------------------------------------
|
||||
! ... reads in and stores the ground state wavefunctions
|
||||
! ... reads in and stores the ground state wavefunctions
|
||||
! ... for use in Lanczos linear response calculation
|
||||
!---------------------------------------------------------------------
|
||||
!
|
||||
|
@ -9,7 +9,7 @@ subroutine lr_read_wf()
|
|||
#include "f_defs.h"
|
||||
!
|
||||
use io_global, only : stdout
|
||||
use klist, only : nks, xk
|
||||
use klist, only : nks, xk
|
||||
use cell_base, only : tpiba2
|
||||
use gvect, only : ngm, g, ecutwfc
|
||||
use io_files, only : nwordwfc, iunwfc, prefix, diropn
|
||||
|
@ -28,17 +28,17 @@ subroutine lr_read_wf()
|
|||
USE realus, ONLY : real_space, fft_orbital_gamma, initialisation_level, &
|
||||
bfft_orbital_gamma, calbec_rs_gamma, add_vuspsir_gamma, &
|
||||
v_loc_psir, s_psir_gamma, real_space_debug
|
||||
USE lr_variables, ONLY : lr_verbosity
|
||||
USE buffers, ONLY : get_buffer
|
||||
USE lr_variables, ONLY : lr_verbosity
|
||||
USE buffers, ONLY : get_buffer
|
||||
USE kinds, ONLY : dp
|
||||
|
||||
|
||||
|
||||
|
||||
!
|
||||
implicit none
|
||||
!
|
||||
!
|
||||
|
||||
|
||||
!
|
||||
! local variables
|
||||
integer :: ik, ibnd, ig, itmp1,itmp2,itmp3
|
||||
|
@ -82,11 +82,11 @@ USE lr_variables, ONLY: check_all_bands_gamma, check_density_gamma,&
|
|||
!
|
||||
if (.not.exst) call errore('lr_read_wfc', TRIM( prefix )//'.wfc'//' not found',1)
|
||||
!
|
||||
if (gamma_only) then
|
||||
if (gamma_only) then
|
||||
WRITE( stdout, '(/5x,"Gamma point algorithm")' )
|
||||
else
|
||||
call errore('lr_read_wfc', 'k-point algorithm is not tested yet',1)
|
||||
WRITE( stdout, '(/5x,"Generalised algorithm !warning")' )
|
||||
WRITE( stdout, '(/5x,"Generalised algorithm !warning")' )
|
||||
endif
|
||||
|
||||
do ik=1,nks
|
||||
|
@ -107,7 +107,7 @@ USE lr_variables, ONLY: check_all_bands_gamma, check_density_gamma,&
|
|||
!
|
||||
!
|
||||
CLOSE( UNIT = iunwfc)
|
||||
!print * , "evc0 ",evc0(1:3,1,1)
|
||||
!print * , "evc0 ",evc0(1:3,1,1)
|
||||
!
|
||||
!
|
||||
! vkb * evc0 and initialization of sevc0
|
||||
|
@ -121,7 +121,7 @@ USE lr_variables, ONLY: check_all_bands_gamma, check_density_gamma,&
|
|||
call init_us_2(npw,igk_k(:,1),xk(1,1),vkb)
|
||||
!
|
||||
if (real_space_debug>0) then
|
||||
!
|
||||
!
|
||||
!
|
||||
!
|
||||
do ibnd=1,nbnd,2
|
||||
|
@ -140,8 +140,8 @@ USE lr_variables, ONLY: check_all_bands_gamma, check_density_gamma,&
|
|||
filename=trim(prefix) // "-rbecp-rs.dump"
|
||||
OPEN(UNIT=47,FILE=filename,STATUS='NEW',ACCESS = 'SEQUENTIAL')
|
||||
write(unit=47,FMT='("#RBECP SIZE :",i6," number of beta fs",i6," bands",i6)') size(becp%r)&
|
||||
,size(becp%r,1),size(becp%r,2)
|
||||
|
||||
,size(becp%r,1),size(becp%r,2)
|
||||
|
||||
do itmp2=1, SIZE(becp%r,2)
|
||||
write(unit=47,FMT='("#Band no",i3)') itmp2
|
||||
do itmp1=1, SIZE(becp%r,1)
|
||||
|
@ -159,7 +159,7 @@ USE lr_variables, ONLY: check_all_bands_gamma, check_density_gamma,&
|
|||
write(unit=48,FMT='(i6,2x,e21.15, 2x, e21.15,2x)') itmp1, DBLE(sevc0(itmp1,itmp2,1)), AIMAG(sevc0(itmp1,itmp2,1))
|
||||
enddo
|
||||
enddo
|
||||
|
||||
|
||||
close(48)
|
||||
endif
|
||||
!print *, becp1-rbecp
|
||||
|
@ -168,19 +168,19 @@ USE lr_variables, ONLY: check_all_bands_gamma, check_density_gamma,&
|
|||
! call s_psi(npwx, npw_k(1), nbnd, evc0(:,:,1), sevc0(:,:,1))
|
||||
else
|
||||
!
|
||||
!call pw_gemm('Y',nkb,nbnd,npw_k(1),vkb,npwx,evc0,npwx,becp1,nkb)
|
||||
!call pw_gemm('Y',nkb,nbnd,npw_k(1),vkb,npwx,evc0,npwx,becp1,nkb)
|
||||
call calbec(npw_k(1),vkb,evc0(:,:,1),becp1)
|
||||
!
|
||||
becp%r=becp1
|
||||
!
|
||||
call s_psi(npwx, npw_k(1), nbnd, evc0(:,:,1), sevc0(:,:,1))
|
||||
! Test case
|
||||
if (test_case_no .eq. 1) then
|
||||
if (test_case_no .eq. 1) then
|
||||
write(stdout,'(/5x,"Test Case 1, dumping Fourier space calculated rbecp and sevc0",1x)')
|
||||
filename=trim(prefix) // "-rbecp.dump"
|
||||
OPEN(UNIT=47,FILE=filename,STATUS='NEW',ACCESS = 'SEQUENTIAL')
|
||||
write(unit=47,FMT='("#RBECP SIZE :",i6," number of beta fs",i6," bands",i6)') size(becp%r)&
|
||||
,size(becp%r,1),size(becp%r,2)
|
||||
,size(becp%r,1),size(becp%r,2)
|
||||
do itmp2=1, SIZE(becp%r,2)
|
||||
write(unit=47,FMT='("#Band no",i3)') itmp2
|
||||
do itmp1=1, SIZE(becp%r,1)
|
||||
|
@ -189,7 +189,7 @@ USE lr_variables, ONLY: check_all_bands_gamma, check_density_gamma,&
|
|||
enddo
|
||||
close(47)
|
||||
filename=trim(prefix) // "-sevc0.dump"
|
||||
OPEN(UNIT=48,FILE=filename,STATUS='NEW',ACCESS = 'SEQUENTIAL')
|
||||
OPEN(UNIT=48,FILE=filename,STATUS='NEW',ACCESS = 'SEQUENTIAL')
|
||||
write(unit=48,FMT='("#SEVC0 SIZE :",i6," NPW ",i6," BANDS ",i6," DIM3",i6)') size(sevc0), &
|
||||
size(sevc0,1), size(sevc0,2), size(sevc0,3)
|
||||
do itmp2=1, SIZE(sevc0,2)
|
||||
|
@ -261,7 +261,7 @@ USE lr_variables, ONLY: check_all_bands_gamma, check_density_gamma,&
|
|||
!
|
||||
endif
|
||||
!
|
||||
call cft3s(revc0(:,ibnd,1),nr1s,nr2s,nr3s,nrx1s,nrx2s,nrx3s,2)
|
||||
CALL invfft ('Wave', revc0(:,ibnd,1), dffts)
|
||||
!
|
||||
end do
|
||||
!
|
||||
|
@ -277,7 +277,7 @@ USE lr_variables, ONLY: check_all_bands_gamma, check_density_gamma,&
|
|||
!
|
||||
enddo
|
||||
!
|
||||
call cft3s(revc0(:,ibnd,ik),nr1s,nr2s,nr3s,nrx1s,nrx2s,nrx3s,2)
|
||||
CALL invfft ('Wave', revc0(:,ibnd,ik), dffts)
|
||||
!
|
||||
end do
|
||||
!
|
||||
|
@ -285,7 +285,7 @@ USE lr_variables, ONLY: check_all_bands_gamma, check_density_gamma,&
|
|||
!
|
||||
end if
|
||||
!
|
||||
!print * , "evc0 ",evc0(1:3,1,1)
|
||||
!print * , "evc0 ",evc0(1:3,1,1)
|
||||
!
|
||||
if (lr_verbosity >10) then
|
||||
call check_all_bands_gamma(evc0(:,:,1),sevc0(:,:,1),nbnd,nbnd)
|
||||
|
@ -322,7 +322,7 @@ USE lr_variables, ONLY: check_all_bands_gamma, check_density_gamma,&
|
|||
real(kind=dp), allocatable :: becp1_all(:,:)
|
||||
complex(kind=dp), allocatable :: becp1_c_all(:,:,:)
|
||||
complex(kind=dp), allocatable :: revc_all(:,:,:)
|
||||
|
||||
|
||||
!First pretend everything is normal
|
||||
nbnd=nbnd_total
|
||||
!
|
||||
|
@ -335,11 +335,11 @@ USE lr_variables, ONLY: check_all_bands_gamma, check_density_gamma,&
|
|||
becp1_all(:,:)=0.0d0
|
||||
else
|
||||
allocate(becp1_c_all(nkb,nbnd,nks))
|
||||
becp1_c_all(:,:,:)=(0.0d0,0.0d0)
|
||||
becp1_c_all(:,:,:)=(0.0d0,0.0d0)
|
||||
endif
|
||||
endif
|
||||
|
||||
|
||||
|
||||
nwordwfc = 2 * nbnd * npwx
|
||||
size_evc=npwx*nbnd_occ(1)*nks
|
||||
!
|
||||
|
@ -347,11 +347,11 @@ USE lr_variables, ONLY: check_all_bands_gamma, check_density_gamma,&
|
|||
!
|
||||
if (.not.exst) call errore('lr_read_wfc', TRIM( prefix )//'.wfc'//' not found',1)
|
||||
!
|
||||
if (gamma_only) then
|
||||
if (gamma_only) then
|
||||
WRITE( stdout, '(/5x,"Gamma point algorithm")' )
|
||||
else
|
||||
call errore('lr_read_wfc', 'k-point algorithm is not tested yet',1)
|
||||
WRITE( stdout, '(/5x,"Generalised algorithm !warning")' )
|
||||
WRITE( stdout, '(/5x,"Generalised algorithm !warning")' )
|
||||
endif
|
||||
|
||||
do ik=1,nks
|
||||
|
@ -372,7 +372,7 @@ USE lr_variables, ONLY: check_all_bands_gamma, check_density_gamma,&
|
|||
!
|
||||
!
|
||||
CLOSE( UNIT = iunwfc)
|
||||
!print * , "evc_all ",evc_all(1:3,1,1)
|
||||
!print * , "evc_all ",evc_all(1:3,1,1)
|
||||
!
|
||||
!
|
||||
! vkb * evc_all and initialization of sevc_all
|
||||
|
@ -386,7 +386,7 @@ USE lr_variables, ONLY: check_all_bands_gamma, check_density_gamma,&
|
|||
call init_us_2(npw,igk_k(:,1),xk(1,1),vkb)
|
||||
!
|
||||
if (real_space_debug>0) then
|
||||
!
|
||||
!
|
||||
!
|
||||
!
|
||||
do ibnd=1,nbnd,2
|
||||
|
@ -463,7 +463,7 @@ USE lr_variables, ONLY: check_all_bands_gamma, check_density_gamma,&
|
|||
!
|
||||
endif
|
||||
!
|
||||
call cft3s(revc_all(:,ibnd,1),nr1s,nr2s,nr3s,nrx1s,nrx2s,nrx3s,2)
|
||||
CALL invfft ('Wave', revc_all(:,ibnd,1), dffts)
|
||||
!
|
||||
end do
|
||||
!
|
||||
|
@ -479,7 +479,7 @@ USE lr_variables, ONLY: check_all_bands_gamma, check_density_gamma,&
|
|||
!
|
||||
enddo
|
||||
!
|
||||
call cft3s(revc_all(:,ibnd,ik),nr1s,nr2s,nr3s,nrx1s,nrx2s,nrx3s,2)
|
||||
CALL invfft ('Wave', revc_all(:,ibnd,ik), dffts)
|
||||
!
|
||||
end do
|
||||
!
|
||||
|
@ -498,7 +498,7 @@ USE lr_variables, ONLY: check_all_bands_gamma, check_density_gamma,&
|
|||
revc0=(0.0d0,0.0d0)
|
||||
revc0(:,:,:)=revc_all(:,1:nbnd,:)
|
||||
if (nkb>0) then
|
||||
if (gamma_only) then
|
||||
if (gamma_only) then
|
||||
becp1(:,:)=becp1_all(:,1:nbnd)
|
||||
becp%r=0.0d0
|
||||
becp%r=becp1
|
||||
|
@ -512,7 +512,7 @@ USE lr_variables, ONLY: check_all_bands_gamma, check_density_gamma,&
|
|||
evc0_virt(:,:,:)=evc_all(:,nbnd+1:nbnd_total,:)
|
||||
!sevc0_virt(:,:,:)=sevc_all(:,nbnd+1:nbnd_total,:)
|
||||
if (nkb>0) then
|
||||
if (gamma_only) then
|
||||
if (gamma_only) then
|
||||
becp1_virt(:,:)=becp1_all(:,nbnd+1:nbnd_total)
|
||||
else
|
||||
becp1_c_virt(:,:,:)=becp1_c_all(:,nbnd+1:nbnd_total,:)
|
||||
|
@ -528,7 +528,7 @@ USE lr_variables, ONLY: check_all_bands_gamma, check_density_gamma,&
|
|||
enddo
|
||||
endif
|
||||
if (nkb>0) then
|
||||
if (gamma_only) then
|
||||
if (gamma_only) then
|
||||
deallocate(becp1_all)
|
||||
else
|
||||
deallocate(becp1_c_all)
|
||||
|
@ -545,4 +545,4 @@ USE lr_variables, ONLY: check_all_bands_gamma, check_density_gamma,&
|
|||
!
|
||||
end subroutine virt_read
|
||||
!-----------------------------------------------------------------------
|
||||
end subroutine lr_read_wf
|
||||
end subroutine lr_read_wf
|
||||
|
|
Loading…
Reference in New Issue