TABLE OF CONTENTS
ABINIT/waveformat [ Functions ]
NAME
waveformat
FUNCTION
This routine is to find the matched pairs of plane waves between two neighbouring k points and load a new pw coefficients array cg_new Was written first by Na Sai (thanks), but unfortunately without any comment ...
COPYRIGHT
Copyright (C) 2000-2018 ABINIT group (NSAI,XG,MKV) 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
cg(2,mpw*nspinor*mband*mkmem*nsppol)= input planewave coefficients, in case mkmem/=0 cg_disk(2,mpw*nspinor*mband,2)= input planewave coefficients, in case mkmem==0 cg_index(mband,nkpt_,nsppol)=index of wavefunction iband,ikpt,isppol in the array cg. dk(3)= step taken to the next k mesh point along the kberry direction (see also isgn) ii=(to be documented) ikpt=index of the first k-point in the reduced Brillouin zone ikpt_=index of the first k-point in the full Brillouin zone isgn=1 if dk(3) is connecting the k-points (ikpt_ and jkpt) =-1 if -dk(3) is connecting the k-points isppol=1 if spin-up, =2 if spin-down jj=(to be documented) jkpt=index of the second k-point in the reduced Brillouin zone jkpt_=index of the second k-point in the full Brillouin zone kg_kpt(:,:,:)= unpacked reduced planewave coordinates with subscript of planewave and k point kpt(3,nkpt)=reduced coordinates of k-point grid that samples the whole BZ kg_jl(3,mpw,2)=(to be documented) maxband/minband= control the minimum and maximum band calculated in the overlap matrix mband=maximum number of bands (dimension of several cg* arrays) mcg=size of wave-functions array (cg) =mpw*nspinor*mband*mkmem*nsppol mcg_disk=size of wave-functions array (cg_disk) =mpw*nspinor*mband mkmem= if 0, the wavefunctions are input in cg_disk, otherwise in cg mpw=maximum number of planewaves (dimension of several cg* arrays) nkpt=number of k points (full Brillouin zone !?!) nkpt_=number of k points (reduced Brillouin zone !?!) npwarr(nkpt)=number of planewaves in basis at this k point nspinor=number of spinorial components of the wavefunctions nsppol=1 for unpolarized, 2 for spin-polarized shift_g_2(nkpt,nkpt)=non-zero if a G vector along kberry is needed to connect k points tr(2)=variable that changes k to -k G to -G $c_g$ to $c_g^*$ when time-reversal symetrie is used
OUTPUT
cg_new(2,mpw,maxband)=planewave coefficients transferred onto the set of planewaves at k
SIDE EFFECTS
PARENTS
uderiv
CHILDREN
SOURCE
68 #if defined HAVE_CONFIG_H 69 #include "config.h" 70 #endif 71 72 #include "abi_common.h" 73 74 75 subroutine waveformat(cg,cg_disk,cg_index,cg_new,dk,ii,ikpt,& 76 & ikpt_,isgn,isppol,jj,jkpt,jkpt_,kg_kpt,kpt,kg_jl,maxband,mband,mcg,mcg_disk,& 77 & minband,mkmem,mpw,nkpt,nkpt_,npwarr,nsppol,nspinor,shift_g_2,tr) 78 79 use m_profiling_abi 80 81 use defs_basis 82 83 !This section has been created automatically by the script Abilint (TD). 84 !Do not modify the following lines by hand. 85 #undef ABI_FUNC 86 #define ABI_FUNC 'waveformat' 87 !End of the abilint section 88 89 implicit none 90 91 !Arguments ------------------------------------ 92 !scalars 93 integer,intent(in) :: ii,ikpt,ikpt_,isgn,isppol,jj,jkpt,jkpt_,maxband,mband,mcg,mcg_disk 94 integer,intent(in) :: minband,mkmem,mpw,nkpt,nkpt_,nspinor,nsppol 95 !arrays 96 integer,intent(in) :: cg_index(mband,nkpt_,nsppol),kg_jl(3,mpw,2) 97 integer,intent(in) :: kg_kpt(3,mpw*nspinor,nkpt_),npwarr(nkpt_) 98 real(dp),intent(in) :: cg(2,mcg) 99 real(dp),intent(in) :: cg_disk(2,mcg_disk,2),dk(3),kpt(3,nkpt),tr(2) 100 real(dp),intent(out) :: cg_new(2,mpw,maxband) 101 logical,intent(in) :: shift_g_2(nkpt,nkpt) 102 103 !Local variables ------------------------- 104 !scalars 105 integer :: cg_index_iband,iband,ipw,jpw,nomatch,npw_k 106 logical :: found_match 107 !arrays 108 integer :: dg(3) 109 110 ! *********************************************************************** 111 112 npw_k=npwarr(ikpt) 113 114 115 nomatch=0 116 117 !If there is no shift of G-vector between ikpt_ and jkpt_ 118 if(shift_g_2(ikpt_,jkpt_) .eqv. .false.) then 119 120 ! DEBUG 121 ! write(111,*)'pair', ikpt_,jkpt_,'noshift' 122 ! ENDDEBUG 123 124 ! If the original wavefunction is contained in cg_disk 125 if(mkmem==0) then 126 127 do ipw=1,npw_k 128 129 found_match = .false. 130 131 do jpw=1,npwarr(jkpt) 132 if (sum(abs(tr(ii)*kg_jl(:,ipw,ii)-tr(jj)*kg_jl(:,jpw,jj)))<3*tol8)then 133 do iband=minband, maxband 134 cg_index_iband=(iband-1)*npwarr(jkpt) 135 cg_new(1:2,ipw,iband)=cg_disk(1:2,jpw+cg_index_iband,jj) 136 end do 137 found_match = .true. 138 exit 139 end if 140 end do 141 142 if (found_match .eqv. .false.) then 143 do iband=minband,maxband 144 cg_new(1:2,ipw,iband)=zero 145 end do 146 nomatch = nomatch + 1 147 end if 148 149 end do 150 151 ! Here, the wavefunctions are contained in cg 152 else 153 154 do ipw=1,npw_k 155 156 found_match = .false. 157 158 do jpw=1,npwarr(jkpt) 159 if (sum(abs(tr(ii)*kg_kpt(:,ipw,ikpt)-tr(jj)*kg_kpt(:,jpw,jkpt)))<3*tol8)then 160 do iband=minband, maxband 161 cg_index_iband=cg_index(iband,jkpt,isppol) 162 cg_new(1:2,ipw,iband)=cg(1:2,jpw+cg_index_iband) 163 end do 164 found_match = .true. 165 exit 166 end if 167 end do 168 169 if (found_match .eqv. .false.) then 170 do iband=minband,maxband 171 cg_new(1:2,ipw,iband)=(0.0_dp,0.0_dp) 172 end do 173 nomatch = nomatch + 1 174 end if 175 end do 176 177 end if 178 179 ! DEBUG 180 ! write(111,*) 'normal pair nomatch=',nomatch 181 ! ENDDEBUG 182 183 ! If there is a G-vector shift between ikpt_ and jkpt_ 184 else 185 186 ! DEBUG 187 ! write(111,*) 'pair',ikpt_,jkpt_,' need shift' 188 ! ENDDEBUG 189 190 dg(:) = -1*nint(tr(jj)*kpt(:,jkpt)-tr(ii)*kpt(:,ikpt)+isgn*dk(:)) 191 192 ! If the original wavefunction is contained in cg_disk 193 if(mkmem==0) then 194 195 do ipw=1,npw_k 196 197 found_match = .false. 198 199 do jpw=1,npwarr(jkpt) 200 if (sum(abs(tr(ii)*kg_jl(:,ipw,ii)-(tr(jj)*kg_jl(:,jpw,jj)-& 201 & dg(:))))<3*tol8)then 202 203 do iband=minband, maxband 204 cg_index_iband=(iband-1)*npwarr(jkpt) 205 cg_new(1:2,ipw,iband)=cg_disk(1:2,jpw+cg_index_iband,jj) 206 end do 207 found_match = .true. 208 exit 209 end if 210 end do 211 212 if (found_match .eqv. .false.) then 213 do iband=minband,maxband 214 cg_new(1:2,ipw,iband)=(0.0_dp,0.0_dp) 215 end do 216 nomatch = nomatch + 1 217 end if 218 end do 219 220 ! Here, the wavefunctions are contained in cg 221 else 222 223 do ipw=1,npw_k 224 225 found_match = .false. 226 227 do jpw=1,npwarr(jkpt) 228 if (sum(abs(tr(ii)*kg_kpt(:,ipw,ikpt)-(tr(jj)*kg_kpt(:,jpw,jkpt)-& 229 & dg(:))))<3*tol8)then 230 do iband=minband, maxband 231 cg_index_iband=cg_index(iband,jkpt,isppol) 232 cg_new(1:2,ipw,iband)=cg(1:2,jpw+cg_index_iband) 233 end do 234 found_match = .true. 235 exit 236 end if 237 end do 238 239 if (found_match .eqv. .false.) then 240 do iband=minband,maxband 241 cg_new(1:2,ipw,iband)=zero 242 end do 243 nomatch = nomatch + 1 244 end if 245 end do 246 247 end if 248 249 ! DEBUG 250 ! write(111,*) 'special pair nomatch=',nomatch 251 ! ENDDEBUG 252 253 end if 254 255 end subroutine waveformat