TABLE OF CONTENTS
ABINIT/symkpt [ Functions ]
NAME
symkpt
FUNCTION
Determines the weights of the k-points for sampling the Brillouin Zone, starting from a first set of weights wtk, and folding it to a new set, by taking into account the symmetries described by symrec, and eventually the time-reversal symmetry. Also compute the number of k points in the reduced set This routine is also used for sampling the q vectors in the Brillouin zone for the computation of thermodynamical properties (from the routine thm9).
COPYRIGHT
Copyright (C) 1999-2018 ABINIT group (XG,LSI) 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
chksymbreak= if 1, will check whether the k point grid is symmetric, and stop if not. gmet(3,3)=reciprocal space metric (bohr**-2). iout=if non-zero, output the new number of kpoints on unit iout kbz(3,nkbz)= k vectors in the BZ. nkbz = number of k-points whose weights are wtk nsym=number of space group symmetries symrec(3,3,nsym)=3x3 matrices of the group symmetries (reciprocal space) timrev: if 1, the time reversal operation has to be taken into account if 0, no time reversal symmetry. wtk(nkbz)=weight assigned to each k point.
OUTPUT
ibz2bz(nkbz)=non-symmetrized indices of the k-points a.k.a. ibz2bz mapping The correspondence beween the iq_ibz point in IBZ and the iq_bz point in the full BZ is obtained via: do ik_ibz=1,nkibz ik_bz = ibz2bz(ik_ibz) end do nkibz = number of k-points in the irreducible set wtk_folded(nkbz)=weight assigned to each k point, taking into account the symmetries
NOTES
The decomposition of the symmetry group in its primitives might speed up the execution. The output variables are stored only in the range 1:nkbz
PARENTS
dfpt_looppert,elphon,ep_setupqpt,get_npert_rbz,getkgrid,harmonic_thermo m_bz_mesh,m_ifc,m_sigmaph
CHILDREN
sort_dp,wrtout
SOURCE
57 #if defined HAVE_CONFIG_H 58 #include "config.h" 59 #endif 60 61 #include "abi_common.h" 62 63 subroutine symkpt(chksymbreak,gmet,ibz2bz,iout,kbz,nkbz,nkibz,nsym,& 64 & symrec,timrev,wtk,wtk_folded) 65 66 use defs_basis 67 use m_profiling_abi 68 use m_errors 69 use m_sort 70 71 !This section has been created automatically by the script Abilint (TD). 72 !Do not modify the following lines by hand. 73 #undef ABI_FUNC 74 #define ABI_FUNC 'symkpt' 75 use interfaces_14_hidewrite 76 !End of the abilint section 77 78 implicit none 79 80 !Arguments ------------------------------- 81 !scalars 82 integer,intent(in) :: chksymbreak,iout,nkbz,nsym,timrev 83 integer,intent(out) :: nkibz 84 !arrays 85 integer,intent(in) :: symrec(3,3,nsym) 86 integer,intent(inout) :: ibz2bz(nkbz) !vz_i 87 real(dp),intent(in) :: gmet(3,3),kbz(3,nkbz),wtk(nkbz) 88 real(dp),intent(out) :: wtk_folded(nkbz) 89 90 !Local variables ------------------------- 91 !scalars 92 integer :: identi,ii,ikpt,ikpt2,ind_ikpt,ind_ikpt2 93 integer :: ikpt_current_length,isym,itim,jj,nkpout,quit,tident 94 real(dp) :: difk,difk1,difk2,difk3,length2trial,reduce,reduce1,reduce2,reduce3 95 character(len=500) :: message 96 !arrays 97 integer,allocatable :: list(:) 98 real(dp) :: gmetkpt(3),ksym(3) 99 real(dp),allocatable :: length2(:) 100 101 ! ********************************************************************* 102 103 if (timrev/=1 .and. timrev/=0) then 104 write(message,'(a,i0)')' timrev should be 0 or 1, while it is equal to ',timrev 105 MSG_BUG(message) 106 end if 107 108 if(nsym/=1)then 109 ! Find the identity symmetry operation 110 do isym=1,nsym 111 tident=1 112 do jj=1,3 113 if(symrec(jj,jj,isym)/=1)tident=0 114 do ii=1,3 115 if( ii/=jj .and. symrec(ii,jj,isym)/=0)tident=0 116 end do 117 end do 118 if(tident==1)then 119 identi=isym 120 !write(message, '(a,i0)' )' symkpt : found identity, with number',identi 121 !call wrtout(std_out,message,'COLL') 122 exit 123 end if 124 end do 125 if(tident==0)then 126 MSG_BUG(' Did not find the identity operation') 127 end if 128 else 129 identi=1 130 end if 131 132 !Initialise the wtk_folded array using the wtk array : 133 do ikpt=1,nkbz 134 wtk_folded(ikpt)=wtk(ikpt) 135 end do 136 137 !Here begins the serious business 138 139 !If there is some possibility for a change (otherwise, wtk_folded is 140 !correctly initialized to give no change) : 141 if(nkbz/=1 .and. (nsym/=1 .or. timrev==1) )then 142 143 ! Store the length of vectors, but take into account umklapp 144 ! processes by selecting the smallest length of all symmetric vectors 145 ABI_ALLOCATE(length2,(nkbz)) 146 147 148 !MG FIXME: 149 ! Here there's a possible problem with the order of symmetries because 150 ! in listkk, time-reversal is the outermost loop. This can create inconsistencies in the symmetry tables. 151 do ikpt=1,nkbz 152 do isym=1,nsym 153 do itim=1,(1-2*timrev),-2 154 ! Get the symmetric of the vector 155 do ii=1,3 156 ksym(ii)=itim*( kbz(1,ikpt)*symrec(ii,1,isym)& 157 & +kbz(2,ikpt)*symrec(ii,2,isym)& 158 & +kbz(3,ikpt)*symrec(ii,3,isym) ) 159 ksym(ii)=ksym(ii)-anint(ksym(ii)+tol8*half) 160 end do 161 gmetkpt(:)=gmet(:,1)*ksym(1)+gmet(:,2)*ksym(2)+gmet(:,3)*ksym(3) 162 length2trial=ksym(1)*gmetkpt(1)+ksym(2)*gmetkpt(2)+ksym(3)*gmetkpt(3) 163 if(isym==1 .and. itim==1)then 164 length2(ikpt)=length2trial 165 else 166 if(length2(ikpt)>length2trial)length2(ikpt)=length2trial 167 end if 168 end do 169 end do 170 end do 171 172 ! Sort the lengths 173 ABI_ALLOCATE(list,(nkbz)) 174 list(:)=(/ (ikpt,ikpt=1,nkbz) /) 175 call sort_dp(nkbz,length2,list,tol14) 176 ! do ikpt=1,nkbz 177 ! write(std_out,*)ikpt,length2(ikpt),list(ikpt),kbz(1:3,list(ikpt)) 178 ! end do 179 180 ! Examine whether the k point grid is symmetric or not 181 if(chksymbreak==1)then 182 ikpt_current_length=1 183 ! Loop on all k points 184 do ikpt=1,nkbz 185 ind_ikpt=list(ikpt) 186 ! Keep track of the current length, to avoid doing needless comparisons 187 if(length2(ikpt)-length2(ikpt_current_length)>tol8)then 188 ikpt_current_length=ikpt 189 end if 190 191 do isym=1,nsym 192 do itim=1,(1-2*timrev),-2 193 if(isym/=identi .or. itim/=1 )then 194 195 ! Get the symmetric of the vector 196 do ii=1,3 197 ksym(ii)=itim*( kbz(1,ind_ikpt)*symrec(ii,1,isym)& 198 & +kbz(2,ind_ikpt)*symrec(ii,2,isym)& 199 & +kbz(3,ind_ikpt)*symrec(ii,3,isym) ) 200 end do 201 202 ! Search over k-points with the same length, to find whether there is a connecting symmetry operation 203 quit=0 204 do ikpt2=ikpt_current_length,nkbz 205 ! The next line skip all ikpt2 vectors, as soon as one becomes larger than length2(ikpt) 206 ! Indeed, one is already supposed to have found a symmetric k point before this happens ... 207 if(length2(ikpt2)-length2(ikpt)>tol8)exit 208 ! Ordered index 209 ind_ikpt2=list(ikpt2) 210 difk1= ksym(1)-kbz(1,ind_ikpt2) 211 reduce1=difk1-anint(difk1) 212 difk2= ksym(2)-kbz(2,ind_ikpt2) 213 reduce2=difk2-anint(difk2) 214 difk3= ksym(3)-kbz(3,ind_ikpt2) 215 reduce3=difk3-anint(difk3) 216 if(abs(reduce1)+abs(reduce2)+abs(reduce3)<tol8)then 217 ! The symmetric was found 218 quit=1 219 exit 220 end if 221 end do 222 if(quit==0)then 223 write(message,'(3a,i4,2a,9i3,2a,i6,1a,3es16.6,6a)' )& 224 & 'Chksymbreak=1 . It has been observed that the k point grid is not symmetric :',ch10,& 225 & 'for the symmetry number ',isym,ch10,& 226 & 'with symrec=',symrec(1:3,1:3,isym),ch10,& 227 & 'the symmetric of the k point number ',ind_ikpt2,' with components', kbz(1:3,ind_ikpt2),ch10,& 228 & 'does not belong to the k point grid.',ch10,& 229 & 'Read the description of the input variable chksymbreak,',ch10,& 230 & 'You might switch it to zero, or change your k point grid to one that is symmetric.' 231 MSG_ERROR(message) 232 end if 233 234 end if ! End condition of non-identity symmetry 235 end do ! itim 236 end do ! isym 237 238 end do ! ikpt 239 end if ! chksymbreak==1 240 241 242 ! Eliminate the k points that are symmetric of another one 243 do ikpt=1,nkbz-1 244 245 ! Ordered index 246 ind_ikpt=list(ikpt) 247 248 ! Not worth to examine a k point that is a symmetric of another, 249 ! which is the case if its weight has been set to 0 by previous folding 250 if(wtk_folded(ind_ikpt)<tol16)cycle 251 252 ! Loop on the remaining k-points 253 do ikpt2=ikpt+1,nkbz 254 255 ! The next line eliminates pairs of vectors that differs by their length. 256 ! Moreover, since the list is ordered according to the length, 257 ! one can skip all other ikpt2 vectors, as soon as one becomes larger than length2(ikpt) 258 if(length2(ikpt2)-length2(ikpt)>tol8)exit 259 260 ! Ordered index 261 ind_ikpt2=list(ikpt2) 262 263 ! If the second vector is already empty, no interest to treat it 264 if(wtk_folded(ind_ikpt2)<tol16)cycle 265 266 quit=0 267 do isym=1,nsym 268 do itim=1,(1-2*timrev),-2 269 if(isym/=identi .or. itim/=1 )then 270 271 ! Get the symmetric of the vector 272 do ii=1,3 273 ksym(ii)=itim*( kbz(1,ind_ikpt)*symrec(ii,1,isym)& 274 & +kbz(2,ind_ikpt)*symrec(ii,2,isym)& 275 & +kbz(3,ind_ikpt)*symrec(ii,3,isym) ) 276 end do 277 278 ! The do-loop was expanded to speed up the execution 279 difk= ksym(1)-kbz(1,ind_ikpt2) 280 reduce=difk-anint(difk) 281 if(abs(reduce)>tol8)cycle 282 difk= ksym(2)-kbz(2,ind_ikpt2) 283 reduce=difk-anint(difk) 284 if(abs(reduce)>tol8)cycle 285 difk= ksym(3)-kbz(3,ind_ikpt2) 286 reduce=difk-anint(difk) 287 if(abs(reduce)>tol8)cycle 288 289 ! Here, have successfully found a symmetrical k-vector 290 ! Assign all the weight of the k-vector to its symmetrical 291 wtk_folded(ind_ikpt)=wtk_folded(ind_ikpt)+wtk_folded(ind_ikpt2) 292 wtk_folded(ind_ikpt2)=0._dp 293 294 ! Go to the next ikpt2 if the symmetric was found 295 quit=1 296 exit 297 298 end if ! End condition of non-identity symmetry 299 end do ! End loop on itim 300 301 if(quit==1)exit 302 end do ! End loop on isym 303 end do ! End secondary loop over k-points 304 end do ! End primary loop over k-points 305 306 ABI_DEALLOCATE(length2) 307 ABI_DEALLOCATE(list) 308 end if ! End check on possibility of change 309 310 !Create the indexing array ibz2bz 311 nkibz=0 312 do ikpt=1,nkbz 313 if(wtk_folded(ikpt)>tol8)then 314 nkibz=nkibz+1 315 ibz2bz(nkibz)=ikpt 316 end if 317 end do 318 319 if(iout/=0)then 320 if(nkbz/=nkibz)then 321 write(message, '(a,a,a,i6,a)' )& 322 & ' symkpt : the number of k-points, thanks to the symmetries,',ch10,& 323 & ' is reduced to',nkibz,' .' 324 call wrtout(iout,message,'COLL') 325 if(iout/=std_out) then 326 call wrtout(std_out,message,'COLL') 327 end if 328 329 nkpout=nkibz 330 !if(nkibz>80)then 331 ! write(message,'(a)' )' greater than 80, so only write 20 of them ' 332 ! call wrtout(std_out,message,'COLL') 333 ! nkpout=20 334 !end if 335 !do ii=1,nkpout 336 ! write(message, '(1x,i2,a2,3es16.8)' ) ii,') ',kbz(1:3,ibz2bz(ii)) 337 ! call wrtout(std_out,message,'COLL') 338 !end do 339 340 ! DEBUG 341 ! write(message, '(a)' )' Here are the new weights :' 342 ! call wrtout(std_out,message,'COLL') 343 ! do ikpt=1,nkbz,6 344 ! write(message, '(6f12.6)' ) wtk_folded(ikpt:min(nkbz,ikpt+5)) 345 ! call wrtout(std_out,message,'COLL') 346 ! end do 347 ! ENDDEBUG 348 else 349 write(message, '(a)' )' symkpt : not enough symmetry to change the number of k points.' 350 call wrtout(iout,message,'COLL') 351 if(iout/=std_out) then 352 call wrtout(std_out,message,'COLL') 353 end if 354 end if 355 end if 356 357 !DEBUG 358 !write(std_out,*)' exit symkpt ' 359 !write(std_out,*)' nkbz,nsym',nkbz,nsym 360 !write(std_out,*)' wtk',wtk 361 !write(std_out,*)' timrev',timrev 362 !if(timrev==0)stop 363 !if(option==1)stop 364 !ENDDEBUG 365 366 end subroutine symkpt