TABLE OF CONTENTS
ABINIT/wfd_vnlpsi [ Functions ]
NAME
wfd_vnlpsi
FUNCTION
Evaluates Vnl |psi> in reciprocal space.
COPYRIGHT
Copyright (C) 2010-2018 ABINIT group (MG) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt . For the initials of contributors, see ~abinit/doc/developers/contributors.txt .
INPUTS
Wfd<wfd_t>=Datatype gathering info on the wavefunctions band=Band index. ik_ibz=K-point index. spin=Spin index npw_k=Number of PW for this k-point (used to dimension output arrays) Cryst<crystal_t>= data type gathering info on symmetries and unit cell GS_hamk<gs_hamiltonian_type)=Information about the Hamiltonian for this k-point. Note that no check is done for the consistency btw ik_ibz and the content of GS_hamk. [Kext]<Kdata_t>=Datatype storing form factors and tables that depend of the K. In the standard mode the routine uses the tables stored in Wfd. Kext is used to calculate Vnl e^{ik'r}|u_k> where k is defined by ik_ibz whereas k' is the k-point that has been used to initialize GS_hamk and Kext. This option is used for the Bloch-state-based interpolation.
OUTPUT
vnl_psi(2,npw_k*Wfd%nspinor)= <G+k|V_nl e^{ik'r}|u_k> opaw_psi(2,npw_k*Wfd%nspinor*Wfd%usepaw) <G+k|1+S|Cnk> e^{ik'r}|u_k>
PARENTS
m_shirley
CHILDREN
load_k_hamiltonian,mkkpg,nonlop,pawcprj_alloc,pawcprj_free,wfd_get_cprj
SOURCE
43 #if defined HAVE_CONFIG_H 44 #include "config.h" 45 #endif 46 47 #include "abi_common.h" 48 49 subroutine wfd_vnlpsi(Wfd,band,ik_ibz,spin,npw_k,Cryst,GS_hamk,vnl_psi,opaw_psi,Kext) 50 51 use defs_basis 52 use defs_datatypes 53 use m_errors 54 use m_profiling_abi 55 56 use m_crystal, only : crystal_t 57 use m_wfd, only : wfd_t, kdata_t, wfd_get_cprj 58 use m_hamiltonian, only : gs_hamiltonian_type, load_k_hamiltonian 59 use m_pawcprj, only : pawcprj_type, pawcprj_alloc, pawcprj_free 60 61 !This section has been created automatically by the script Abilint (TD). 62 !Do not modify the following lines by hand. 63 #undef ABI_FUNC 64 #define ABI_FUNC 'wfd_vnlpsi' 65 use interfaces_66_nonlocal 66 !End of the abilint section 67 68 implicit none 69 70 !Arguments ------------------------------------ 71 !scalars 72 integer,intent(in) :: band,ik_ibz,spin,npw_k 73 type(crystal_t),intent(in) :: Cryst 74 type(gs_hamiltonian_type),intent(inout) :: GS_hamk 75 type(Kdata_t),optional,target,intent(in) :: Kext 76 type(wfd_t),target,intent(inout) :: Wfd 77 !arrays 78 real(dp),intent(out) :: vnl_psi(2,npw_k*Wfd%nspinor) ! <G|Vnl|Cnk> 79 real(dp),intent(out) :: opaw_psi(2,npw_k*Wfd%nspinor*Wfd%usepaw) ! <G|1+S|Cnk> 80 81 !Local variables ------------------------------ 82 !scalars 83 integer,parameter :: idir0=0,ider0=0,nnlout0=0,tim_nonlop0=0 84 integer :: choice,cpopt,cp_dim,istwf_k,natom,nkpg,nspinor 85 integer :: paw_opt,signs 86 logical :: use_kptext,ltest 87 character(len=500) :: msg 88 !arrays 89 integer,ABI_CONTIGUOUS pointer :: kg_k(:,:) 90 real(dp) :: kpoint(3),dum_enlout(0),dummy_lambda(1) 91 real(dp),ABI_CONTIGUOUS pointer :: ffnl(:,:,:,:),ph3d(:,:,:) 92 real(dp),allocatable :: kpg_k(:,:),vectin(:,:) 93 type(pawcprj_type),allocatable :: Cprj(:,:) 94 95 !************************************************************************ 96 97 DBG_ENTER("COLL") 98 99 ABI_CHECK(Wfd%nspinor==1,"nspinor==2 not coded") 100 101 use_kptext = PRESENT(Kext) 102 103 natom = Cryst%natom 104 nspinor = Wfd%nspinor 105 106 !npw_k = Wfd%Kdata(ik_ibz)%npw 107 kg_k => Wfd%Kdata(ik_ibz)%kg_k 108 109 if (.not.use_kptext) then 110 istwf_k = Wfd%Kdata(ik_ibz)%istwfk 111 kpoint = Wfd%kibz(:,ik_ibz) 112 ffnl => Wfd%Kdata(ik_ibz)%fnl_dir0der0 113 ph3d => Wfd%Kdata(ik_ibz)%ph3d 114 else 115 istwf_k = Kext%istwfk 116 kpoint = gs_hamk%kpt_k 117 ffnl => Kext%fnl_dir0der0 118 ph3d => Kext%ph3d 119 ltest = (istwf_k==1 .and. Wfd%istwfk(ik_ibz)==1) 120 if (.not.ltest) then 121 write(msg,'(2a,2(a,i0))')& 122 & "istwfk and Wfd%istwfk(ik_ibz) must be 1 when Kext is present, however, ",ch10,& 123 & "istwfk= ",istwf_k," and Wfd%istwfk= ",Wfd%istwfk(ik_ibz) 124 MSG_BUG(msg) 125 end if 126 ltest = (ik_ibz==1 .and. ALL(ABS(Wfd%kibz(:,ik_ibz))<tol6)) 127 if (.not.ltest) then 128 write(msg,'(3a,i0,a,3(f8.3,1x))')& 129 & "ik_ibz must be 1 and kibz should be 0 when Kext is present, however, ",ch10,& 130 & "ik_ibz= ",ik_ibz," and kpoint= ",Wfd%kibz(:,ik_ibz) 131 MSG_BUG(msg) 132 end if 133 end if 134 ! 135 ! Input wavefunction coefficients <G|Cnk> 136 ABI_MALLOC(vectin,(2,npw_k*nspinor)) 137 vectin(1,:) = DBLE (Wfd%Wave(band,ik_ibz,spin)%ug) 138 vectin(2,:) = AIMAG(Wfd%Wave(band,ik_ibz,spin)%ug) 139 140 signs = 2 ! => apply the non-local operator to a function in G-space. 141 choice = 1 ! => <G|V_nonlocal|vectin>. 142 cpopt =-1 143 paw_opt= 0 144 if (Wfd%usepaw==1) then 145 paw_opt=4 ! both PAW nonlocal part of H (Dij) and overlap matrix (Sij) 146 cpopt=3 ! <p_lmn|in> are already in memory 147 if (use_kptext) cpopt=-1 ! Since we have to change the k-point <p_lmn|in> are recomputed in nonlocal and not saved 148 end if 149 150 cp_dim = ((cpopt+5)/5) 151 ABI_DT_MALLOC(Cprj,(natom,nspinor*cp_dim)) 152 153 if (cp_dim>0) then 154 call pawcprj_alloc(Cprj,0,Wfd%nlmn_sort) 155 call wfd_get_cprj(Wfd,band,ik_ibz,spin,Cryst,Cprj,sorted=.TRUE.) 156 end if 157 ! 158 ! Compute (k+G) vectors 159 ! Not needed here, however they might be stored in Kdata_t to avoid the calculation inside the loop over bands 160 !nkpg=3*GS_hamk%nloalg(3) 161 nkpg=0 162 ABI_MALLOC(kpg_k,(npw_k,nkpg)) 163 if (nkpg>0) then 164 call mkkpg(kg_k,kpg_k,kpoint,nkpg,npw_k) 165 endif 166 167 call load_k_hamiltonian(GS_hamk,kpt_k=kpoint,npw_k=npw_k,istwf_k=istwf_k,& 168 & kg_k=kg_k,kpg_k=kpg_k,ffnl_k=ffnl,ph3d_k=ph3d) 169 170 call nonlop(choice,cpopt,Cprj,dum_enlout,GS_hamk,idir0,dummy_lambda,Wfd%MPI_enreg,1,nnlout0,& 171 & paw_opt,signs,opaw_psi,tim_nonlop0,vectin,vnl_psi) 172 173 ABI_FREE(kpg_k) 174 ABI_FREE(vectin) 175 176 if (cp_dim>0) then 177 call pawcprj_free(Cprj) 178 end if 179 ABI_DT_FREE(Cprj) 180 181 DBG_EXIT("COLL") 182 183 end subroutine wfd_vnlpsi