TABLE OF CONTENTS
ABINIT/m_fold2block [ Modules ]
NAME
m_fold2block
FUNCTION
This module contains basic tools to operate on vectors expressed in reduced coordinates.
COPYRIGHT
Copyright (C) 2008-2022 ABINIT group (AB) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
SOURCE
16 #if defined HAVE_CONFIG_H 17 #include "config.h" 18 #endif 19 20 #include "abi_common.h" 21 22 MODULE m_fold2block 23 24 use defs_basis 25 use m_abicore 26 use m_errors 27 28 implicit none 29 30 public 31 !private 32 33 CONTAINS !===========================================================
m_fold2block/getargs [ Functions ]
[ Top ] [ m_fold2block ] [ Functions ]
NAME
getargs
FUNCTION
Part of fold2Bloch that interprets the command line arguments and checks for their corectness.
INPUTS
fname: Name of input file folds: Number of folds in x, y, and z directions.
OUTPUT
folds: Number of folds in x, y, and z directions. Returns the folds array filled with each element.
SOURCE
243 subroutine getargs(folds, fname) 244 245 implicit none 246 247 !Arguments ------------------------------------ 248 !scalars 249 character(len=fnlen), intent(inout) :: fname 250 !arrays 251 integer, intent(inout) :: folds(3) 252 253 !Local variables------------------------------- 254 !scalars 255 integer :: num_args, argcount, ii, ios 256 character(len=fnlen) :: argfolds 257 logical :: dir 258 !arrays 259 character(len=fnlen), allocatable :: args(:) 260 261 ! ************************************************************************* 262 263 !get number of args 264 num_args=command_argument_count() 265 ABI_MALLOC(args,(num_args)) 266 267 !Check for errors in arguments 268 if (num_args>2) then 269 write(std_out,*) " Too many arguments." 270 write(std_out,*) " Usage: $fold2Bloch file_WFK x:y:z (folds)" 271 ABI_ERROR("Aborting now") 272 elseif (num_args<2) then 273 if (num_args==1) then 274 call get_command_argument(1,args(1)) 275 if ((args(1)=="-h").or.(args(1)=="--help")) then 276 write(std_out,*) " Usage: fold2Bloch file_WFK x:y:z (folds)" 277 write(std_out,*) " file_WFK is the WFK file name (ex. Ex_WFK)" 278 write(std_out,*) " x:y:z integers, greater than 0 that represent a multiplicity" 279 write(std_out,*) " in the corresponding directions used when constructing the supercell." 280 ABI_ERROR("Aborting now") 281 else 282 write(std_out,*) " Not all arguments are present." 283 write(std_out,*) " Make sure that file name and number of folds are indicated." 284 write(std_out,*) " Usage: $fold2Bloch file_WFK x:y:z (folds)" 285 ABI_ERROR("Aborting now") 286 end if 287 else 288 write(std_out,*) " Not all arguments are present." 289 write(std_out,*) " Make sure that file name and number of folds are indicated." 290 write(std_out,*) " Usage: $fold2Bloch file_WFK x:y:z (folds)" 291 ABI_ERROR("Aborting now") 292 end if 293 else 294 do argcount=1, num_args 295 call get_command_argument(argcount,args(argcount)) 296 if (argcount==1) then 297 read(args(argcount), "(a)") fname !read first argument 298 else 299 read(args(argcount), "(a)") argfolds !read second argument 300 end if 301 end do 302 end if 303 304 ! Does intput file exist? 305 inquire(file=fname, exist=dir) 306 if (.not.(dir)) then 307 write(std_out,*) " Case file not found: ", trim(fname) 308 write(std_out,*) " Usage: $fold2Bloch file_WFK x:y:z (folds)" 309 ABI_ERROR("Aborting now") 310 end if 311 312 !Was the number of folds entered in correct format? 313 ii=0 314 ii=INDEX(argfolds, ':') !look for first ":" to differentiate between axis 315 if (ii==0) then 316 write(std_out,*) " Unknown number of folds. See below or type:" 317 write(std_out,*) " fold2Bloch <-h> or fold2Bloch <--help> for more information." 318 write(std_out,*) ' Usage: $fold2Bloch file_WFK x:y:z (folds)' 319 ABI_ERROR("Aborting now") 320 end if 321 read (argfolds(1:ii-1), *, iostat=ios) folds(1) !read X folds 322 if ((ios/=0).or.(folds(1)<=0)) then 323 ABI_ERROR('Number of folds has to be a positive integer greater than 0') 324 end if 325 argfolds=argfolds(ii+1:) !Start argfolds from the first ":" 326 ii=0 327 ii=INDEX(argfolds, ':') !look for second ":" 328 if (ii==0) then 329 write(std_out,*) ' Unknown number of folds. See below or type:' 330 write(std_out,*) " fold2Bloch <-h> or fold2Bloch <--help> for more information." 331 write(std_out,*) ' Usage: $fold2Bloch file_WFK x:y:z (folds)' 332 ABI_ERROR("Aborting now") 333 end if 334 read (argfolds(1:ii-1),*, iostat=ios) folds(2) !read Y folds 335 if ((ios/=0).or.(folds(2)<=0)) then 336 ABI_ERROR('Number of folds has to be a positive integer greater than 0') 337 end if 338 read(argfolds(ii+1:),*, iostat=ios) Folds(3) !read Z folds 339 if ((ios/=0).or.(folds(3)<=0)) then 340 ABI_ERROR('Number of folds has to be a positive integer greater than 0') 341 end if 342 343 ABI_FREE(args) 344 345 end subroutine getargs
m_fold2block/NewK [ Functions ]
[ Top ] [ m_fold2block ] [ Functions ]
NAME
NewK
FUNCTION
Part of fold2Bloch that determines the unfolded positions of each K point.
INPUTS
FX, FY, FZ= Number of folds in each dimention XX, YY, ZZ= Original K point coordinates
OUTPUT
nkval= Array of determined coordinates of unfolded K point locaations. Depends on the number of folds the WFK file was structured with.
SOURCE
160 subroutine NewK(XX, YY, ZZ, FX, FY, FZ, NKVal) 161 162 implicit none 163 164 !Arguments ------------------------------------ 165 !scalars 166 integer, intent(in) :: FX, FY, FZ 167 real(dp), intent(in) :: XX, YY, ZZ 168 real(dp), intent(inout) :: NKVal (3,FX*FY*FZ) 169 !arrays 170 171 !Local variables------------------------------- 172 !scalars 173 real(dp) :: field_x, field_y, field_z 174 real(dp) :: TX, TY, TZ 175 integer :: loop, ii, jj, kk, size 176 !arrays 177 178 ! ************************************************************************* 179 180 size=FX*FY*FZ 181 182 !Brillouin zone range 183 field_x=0.5*FX 184 field_y=0.5*FY 185 field_z=0.5*FZ 186 loop=1 187 188 !Determine new K point states 189 do ii=0, (FX-1) 190 if ((XX+ii)>(field_x)) then !if falls outside of field, loop it around 191 TX=XX-(field_x*2)+ii 192 else 193 TX=XX+ii 194 end if 195 do jj=0, (FY-1) 196 if ((YY+jj)>(field_y)) then 197 TY=YY-(field_y*2)+jj 198 else 199 TY=YY+jj 200 end if 201 do kk=0,(FZ-1) 202 if ((ZZ+kk)>(field_z)) then 203 TZ=ZZ-(field_z*2)+kk 204 else 205 TZ=ZZ+kk 206 end if 207 NKVal(1,loop)=TX !Assign to an output array 208 NKVal(2,loop)=TY 209 NKVal(3,loop)=TZ 210 loop=loop+1 211 end do 212 end do 213 end do 214 215 !reduce the values to fit in Brillouin zone 216 do ii=1, FX*FY*FZ 217 nkval(1,ii)=nkval(1,ii)/FX 218 nkval(2,ii)=nkval(2,ii)/FY 219 nkval(3,ii)=nkval(3,ii)/FZ 220 end do 221 222 END SUBROUTINE NewK
m_fold2block/progress [ Functions ]
[ Top ] [ m_fold2block ] [ Functions ]
NAME
progress
FUNCTION
Part of fold2Bloch that unfolds the wavefunction and calculates the weight of each band.
INPUTS
ikpt: K point index nkpt: Total number of K points in the function kpt: Current K point being procesed.
OUTPUT
Write to screen K point being processed and percent complete.
SOURCE
366 subroutine progress(ikpt, nkpt, kpt) 367 368 implicit none 369 370 !Arguments ------------------------------------ 371 !scalars 372 integer, intent(in) :: ikpt, nkpt 373 !arrays 374 real(dp), intent(in) :: kpt(3) 375 376 ! ************************************************************************* 377 378 write(std_out, '(a,i10, a1)', advance='no') ''//achar(27)//'[97m',100*ikpt/nkpt,'%'//achar(27)//'[0m' 379 write(std_out, '(a25,3f12.6,a)') ''//achar(27)//'[32m Processing K point: ', kpt,''//achar(27)//'[0m' 380 381 end subroutine progress
m_fold2block/SortC [ Functions ]
[ Top ] [ m_fold2block ] [ Functions ]
NAME
SortC
FUNCTION
Part of fold2Bloch that unfolds the wavefunction and calculates the weight of each band.
INPUTS
FX, FY, FZ= Number of folds in each dimention vector= Array of energy coefficient position vectors coefc= Array of energy coefficients of in a certain band NV= Number of vectors/coefficients
OUTPUT
wegihts= Calculated weights of a band in an unfolded state. Depends on the number of folds the WFK file was structured with.
SOURCE
59 subroutine SortC(FX, FY, FZ, Vector, CoefC, NV, Weights) 60 61 implicit none 62 63 !Arguments ------------------------------------ 64 !scalars 65 integer, intent(in) :: FX, FY, FZ, NV 66 !arrays 67 integer, intent(in) :: Vector(3,NV) 68 real(dp), intent(in) :: CoefC(2,NV) 69 real(dp), intent(inout) :: Weights(FX*FY*FZ) 70 71 !Local variables------------------------------- 72 !scalars 73 real :: sumtot 74 integer :: remainder_x, remainder_y, remainder_z, jj, kk, ll, pp, el 75 !arrays 76 real(dp), allocatable :: Sums(:) 77 real(dp), allocatable :: coefsqr(:),TGroupC(:,:,:,:) 78 integer, allocatable :: counter(:,:,:) 79 80 ! ************************************************************************* 81 82 ABI_MALLOC(TGroupC,(FX, FY, FZ, NV)) 83 ABI_MALLOC(Sums,(FX*FY*FZ)) 84 ABI_MALLOC(counter,(FX,FY,FZ)) 85 ABI_MALLOC(coefsqr,(NV)) 86 87 !Convert to sum of squares 88 do jj=1, NV 89 coefsqr(jj)=coefc(1,jj)*coefc(1,jj)+coefc(2,jj)*coefc(2,jj) 90 end do 91 92 !Initiates the counter and TGroupC elements at 0 93 do jj=1,FX 94 do kk=1,FY 95 do ll=1,FZ 96 counter(jj,kk,ll)=0 97 do pp=1,NV 98 TGroupC(jj,kk,ll,pp)=0.0 99 end do 100 end do 101 end do 102 end do 103 104 !Sorts the energy coeeficients 105 do jj=1, (NV) 106 remainder_x=MODULO(Vector(1,jj), FX) 107 remainder_y=MODULO(Vector(2,jj), FY) 108 remainder_z=MODULO(Vector(3,jj), FZ) 109 counter(remainder_x+1, remainder_y+1, remainder_z+1)=counter(remainder_x+1, remainder_y+1, remainder_z+1)+1 110 TGroupC(remainder_x+1, remainder_y+1, remainder_z+1, counter(remainder_x+1, remainder_y+1, remainder_z+1))=Coefsqr(jj) 111 end do 112 113 !Sums the squares of all coefficients per group 114 el=1 115 do jj=1, FX 116 do kk=1, FY 117 do ll=1, FZ 118 if (counter(jj, kk, ll)>0) then 119 Sums(el)=SUM(TGroupC(jj, kk, ll,1:counter(jj, kk, ll))) 120 el=el+1 121 else 122 Sums(el)=0.0 123 el=el+1 124 end if 125 end do 126 end do 127 end do 128 129 !Assign weights to an array 130 sumtot=SUM(Sums) 131 do jj=1, (FX*FY*FZ) 132 Weights(jj)=Sums(jj)/sumtot 133 end do 134 ABI_FREE(TGroupC) 135 ABI_FREE(Sums) 136 ABI_FREE(counter) 137 ABI_FREE(coefsqr) 138 139 end subroutine SortC