TABLE OF CONTENTS
- ABINIT/m_lgroup
- m_lgroup/lgroup_find_ibzimage
- m_lgroup/lgroup_findq_ibzk
- m_lgroup/lgroup_print
- m_sigmaph/lgroup_free
- m_sigmaph/lgroup_new
- m_sigmaph/lgroup_t
ABINIT/m_lgroup [ Modules ]
NAME
m_lgroup
FUNCTION
The little group of a q-point is defined as the subgroup of rotations that preserves q, modulo a G0 vector (also called umklapp vector). Namely: Sq = q + G0 where S is an operation in reciprocal space (symrec) If time reversal symmetry can be used, it is possible to enlarge the little group by including the operations such as: -Sq = q + G0 The operations of little group define an irriducible wedge, IBZ(q), that is usually larger than the irredubile zone defined by the point group of the crystal. The two zones coincide when q=0
COPYRIGHT
Copyright (C) 2008-2024 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 .
SOURCE
28 #if defined HAVE_CONFIG_H 29 #include "config.h" 30 #endif 31 32 #include "abi_common.h" 33 34 module m_lgroup 35 36 use defs_basis 37 use m_errors 38 use m_abicore 39 use m_crystal 40 use m_copy 41 use m_symkpt 42 use m_sort 43 use m_xmpi 44 45 use m_fstrings, only : ftoa, ktoa, sjoin, ltoa 46 use m_numeric_tools, only : wrap2_pmhalf 47 use m_geometry, only : normv 48 use m_kpts, only : listkk 49 use m_symtk, only : sg_multable, littlegroup_q 50 51 implicit none 52 53 private
m_lgroup/lgroup_find_ibzimage [ Functions ]
[ Top ] [ m_lgroup ] [ Functions ]
NAME
lgroup_find_ibzimage
FUNCTION
Find the symmetrical image in the IBZ(k) of a qpoint in the BZ Returns -1 if not found.
INPUTS
OUTPUT
SOURCE
379 integer function lgroup_find_ibzimage(self, qpt) result(iq_ibz) 380 381 !Arguments ------------------------------------ 382 class(lgroup_t),intent(in) :: self 383 real(dp),intent(in) :: qpt(3) 384 385 !Local variables------------------------------- 386 !scalars 387 integer, parameter :: timrev0 = 0 388 real(dp) :: dksqmax 389 !arrays 390 integer :: indkk(6) 391 ! ************************************************************************* 392 393 ! Note use_symrec and timrev0 394 call listkk(dksqmax, self%gmet, indkk, self%ibz, qpt, self%nibz, 1, self%nsym_lg, & 395 1, self%symafm_lg, self%symrec_lg, timrev0, xmpi_comm_self, use_symrec=.True.) 396 397 iq_ibz = indkk(1) 398 if (dksqmax > tol12) iq_ibz = -1 399 400 end function lgroup_find_ibzimage
m_lgroup/lgroup_findq_ibzk [ Functions ]
[ Top ] [ m_lgroup ] [ Functions ]
NAME
lgroup_findq_ibzk
FUNCTION
Find the index of the q-point in the IBZ(k). Umklapp vectors are not allowed. Returns -1 if not found.
INPUTS
qpt(3)=q-point in reduced coordinates. [qtol]=Optional tolerance for q-point comparison. For each reduced direction the absolute difference between the coordinates must be less that qtol
SOURCE
336 integer pure function lgroup_findq_ibzk(self, qpt, qtol) result(iqpt) 337 338 !Arguments ------------------------------------ 339 !scalars 340 real(dp),optional,intent(in) :: qtol 341 class(lgroup_t),intent(in) :: self 342 !arrays 343 real(dp),intent(in) :: qpt(3) 344 345 !Local variables------------------------------- 346 integer :: iq 347 real(dp) :: my_qtol 348 349 ! ************************************************************************* 350 351 my_qtol = tol6; if (present(qtol)) my_qtol = qtol 352 353 iqpt = -1 354 do iq=1,self%nibz 355 if (all(abs(self%ibz(:, iq) - qpt) < my_qtol)) then 356 iqpt = iq; exit 357 end if 358 end do 359 360 end function lgroup_findq_ibzk
m_lgroup/lgroup_print [ Functions ]
[ Top ] [ m_lgroup ] [ Functions ]
NAME
lgroup_print
FUNCTION
Print the object
INPUTS
[title]=String to be printed as header for additional info. [unit]=Unit number for output [prtvol]=Verbosity level
OUTPUT
Only printing
SOURCE
422 subroutine lgroup_print(self, title, unit, prtvol) 423 424 !Arguments ------------------------------------ 425 integer,optional,intent(in) :: unit, prtvol 426 character(len=*),optional,intent(in) :: title 427 class(lgroup_t),intent(in) :: self 428 429 !Local variables------------------------------- 430 !scalars 431 integer :: my_prtvol, my_unt, ik, ii 432 character(len=500) :: msg 433 ! ************************************************************************* 434 435 my_unt = std_out; if (present(unit)) my_unt = unit 436 my_prtvol = 0; if (present(prtvol)) my_prtvol = prtvol 437 438 msg = ' ==== Info on the <lgroup_t> object ==== ' 439 if (present(title)) msg = ' ==== '//trim(adjustl(title))//' ==== ' 440 call wrtout(my_unt, msg) 441 442 write(msg, '(3a, 2(a, i0, a))') & 443 ' Little group point: ................... ', trim(ktoa(self%point)), ch10, & 444 ' Number of points in IBZ(p) ............ ', self%nibz, ch10, & 445 ' Time-reversal flag (0: No, 1: Yes) .... ', self%input_timrev, ch10 446 call wrtout(my_unt, msg) 447 448 if (my_prtvol > 1) then 449 do ii=1,self%nsym_lg 450 call wrtout(std_out, sjoin("lgsym2glob:", ltoa(self%lgsym2glob(:, ii)))) 451 end do 452 do ik=1,self%nibz 453 call wrtout(my_unt, sjoin(ktoa(self%ibz(:,ik)), ftoa(self%weights(ik)))) 454 end do 455 end if 456 457 end subroutine lgroup_print
m_sigmaph/lgroup_free [ Functions ]
[ Top ] [ m_sigmaph ] [ Functions ]
NAME
lgroup_free
FUNCTION
Free memory
SOURCE
469 subroutine lgroup_free(self) 470 471 !Arguments ------------------------------------ 472 class(lgroup_t),intent(inout) :: self 473 474 ! ************************************************************************* 475 476 ! integer 477 ABI_SFREE(self%symrec_lg) 478 ABI_SFREE(self%symafm_lg) 479 ABI_SFREE(self%bz2ibz_smap) 480 ABI_SFREE(self%lgsym2glob) 481 ABI_SFREE(self%symtab) 482 ! real 483 ABI_SFREE(self%ibz) 484 ABI_SFREE(self%weights) 485 486 end subroutine lgroup_free
m_sigmaph/lgroup_new [ Functions ]
[ Top ] [ m_sigmaph ] [ Functions ]
NAME
lgroup_new
FUNCTION
Build the little group of the k-point. Return IBZ(k) points packed in shells. to facilitate optimization of loops.
INPUTS
cryst(crystal_t)=Crystalline structure kpoint(3)=External k-point defining the little-group timrev=1 if time-reversal symmetry can be used, 0 otherwise. nkbz=Number of k-points in the BZ. kbz(3,nkbz)=K-points in the BZ. nkibz=Number of k-points in the IBZ kibz(3,nkibz)=Irreducible zone. comm= MPI communicator. sord=Defines how to order the points in %ibz. ">" for increasing norm. "<" decreasing. Default: ">"
SOURCE
168 type(lgroup_t) function lgroup_new(cryst, kpoint, timrev, nkbz, kbz, nkibz, kibz, comm, sord) result(new) 169 170 !Arguments ------------------------------------ 171 !scalars 172 integer,intent(in) :: timrev,nkibz,nkbz,comm 173 type(crystal_t),intent(in) :: cryst 174 character(len=1),optional,intent(in) :: sord 175 !arrays 176 real(dp),intent(in) :: kpoint(3),kbz(3,nkbz),kibz(3,nkibz) 177 178 !Local variables ------------------------------ 179 !scalars 180 integer,parameter :: iout0=0, my_timrev0 = 0, chksymbreak0 = 0, debug = 0 181 integer :: otimrev_k,ierr,itim,isym,ik_ibz,ik_bz,ksign,isym_lgk 182 !arrays 183 integer :: symrec_lg(3,3,2*cryst%nsym), symafm_lg(2*cryst%nsym), lgsym2glob(2, 2*cryst%nsym) 184 integer,allocatable :: ibz2bz(:), iperm(:), inv_iperm(:) 185 real(dp) :: kred(3), shift(3) 186 real(dp),allocatable :: wtk_folded(:), kord(:,:) 187 188 ! ************************************************************************* 189 190 ! TODO: Add option to exclude umklapp/time-reversal symmetry and kptopt 191 new%point = kpoint 192 new%input_timrev = timrev 193 new%gmet = cryst%gmet 194 new%nbz = nkbz 195 196 ! Determines the symmetry operations by which the k-point is preserved, 197 ABI_MALLOC(new%symtab, (4, 2, cryst%nsym)) 198 call littlegroup_q(cryst%nsym, kpoint, new%symtab, cryst%symrec, cryst%symafm, otimrev_k, prtvol=0) 199 200 new%nsym_lg = 0 201 do itim=1,new%input_timrev + 1 202 do isym=1,cryst%nsym 203 if (cryst%symafm(isym) == -1) cycle 204 if (new%symtab(4, itim, isym) /= 1) cycle ! not \pm Sq = q+g0 205 new%nsym_lg = new%nsym_lg + 1 206 symrec_lg(:,:,new%nsym_lg) = cryst%symrec(:,:,isym) * (-2* itim + 3) 207 symafm_lg(new%nsym_lg) = cryst%symafm(isym) 208 lgsym2glob(:, new%nsym_lg) = [isym, itim] 209 end do 210 end do 211 212 call alloc_copy(symrec_lg(:, :, 1:new%nsym_lg), new%symrec_lg) 213 call alloc_copy(symafm_lg(1:new%nsym_lg), new%symafm_lg) 214 call alloc_copy(lgsym2glob(:, 1:new%nsym_lg), new%lgsym2glob) 215 216 ! Check group closure. 217 if (debug /= 0) then 218 call sg_multable(new%nsym_lg, symafm_lg, symrec_lg, ierr) 219 ABI_CHECK(ierr == 0, "Error in group closure") 220 end if 221 222 ! Find the irreducible zone with the little group operations. 223 ! Do not use time-reversal since it has been manually introduced previously 224 ABI_MALLOC(ibz2bz, (nkbz)) 225 ABI_MALLOC(new%bz2ibz_smap, (6, nkbz)) 226 ! IBZ2BZ ? 227 228 ! TODO: In principle here we would like to have a set that contains the initial IBZ. 229 call symkpt_new(chksymbreak0, cryst%gmet, ibz2bz, iout0, kbz, nkbz, new%nibz,& 230 new%nsym_lg, new%symrec_lg, my_timrev0, new%bz2ibz_smap, comm) 231 232 ABI_MALLOC(new%ibz, (3, new%nibz)) 233 ABI_CALLOC(new%weights, (new%nibz)) 234 235 do ik_bz=1,nkbz 236 ik_ibz = new%bz2ibz_smap(1, ik_bz) 237 isym_lgk = new%bz2ibz_smap(2, ik_bz) 238 new%bz2ibz_smap(2, ik_bz) = lgsym2glob(1, isym_lgk) 239 new%bz2ibz_smap(3, ik_bz) = lgsym2glob(2, isym_lgk) 240 new%weights(ik_ibz) = new%weights(ik_ibz) + 1 241 end do 242 new%weights(:) = new%weights(:) / nkbz 243 244 do ik_ibz=1,new%nibz 245 ik_bz = ibz2bz(ik_ibz) 246 new%ibz(:, ik_ibz) = kbz(:, ik_bz) 247 end do 248 249 ! TODO: Activate this part so that we can cache the q-point in the IBZ. 250 ! Results are ok but this change is postponed because it leads to an increase 251 ! in the walltime spent in listkk likely because of the different order. 252 253 ! Need to repack the IBZ points and rearrange the other arrays dimensioned with nibz. 254 ! In principle, the best approach would be to pack in stars using crystal%symrec. 255 ! For the time being we pack in shells (much easier). Use wtk_folded as workspace to store the norm. 256 257 !ksign = 0 258 ksign = +1 259 if (present(sord)) then 260 if (sord == "<") ksign = -1 261 if (sord == ">") ksign = +1 262 !if (sord == "r") ksign = 0 263 end if 264 265 if (ksign /= 0) then 266 ABI_MALLOC(wtk_folded, (new%nibz)) 267 do ik_ibz=1,new%nibz 268 call wrap2_pmhalf(new%ibz(:, ik_ibz), kred, shift) 269 wtk_folded(ik_ibz) = ksign * normv(kred, cryst%gmet, "G") 270 end do 271 272 ABI_MALLOC(iperm, (new%nibz)) 273 iperm = [(ik_ibz, ik_ibz=1, new%nibz)] 274 call sort_dp(new%nibz, wtk_folded, iperm, tol12) 275 !iperm = [(ik_ibz, ik_ibz=1, new%nibz)] 276 277 ! Trasfer data. 278 ABI_MALLOC(kord, (3, new%nibz)) 279 do ik_ibz=1,new%nibz 280 kord(:, ik_ibz) = new%ibz(:, iperm(ik_ibz)) 281 wtk_folded(ik_ibz) = new%weights(iperm(ik_ibz)) 282 end do 283 new%ibz = kord(:, 1:new%nibz) 284 new%weights = wtk_folded(1:new%nibz) 285 286 ! Rearrange bz2ibz_smap as well --> need the inverse of iperm. 287 ABI_MALLOC(inv_iperm, (new%nibz)) 288 do ik_ibz=1,new%nibz 289 inv_iperm(iperm(ik_ibz)) = ik_ibz 290 end do 291 292 do ik_bz=1,new%nbz 293 ik_ibz = new%bz2ibz_smap(1, ik_bz) 294 new%bz2ibz_smap(1, ik_bz) = inv_iperm(ik_ibz) 295 end do 296 297 ABI_FREE(inv_iperm) 298 ABI_FREE(kord) 299 ABI_FREE(iperm) 300 ABI_FREE(wtk_folded) 301 end if 302 303 ABI_FREE(ibz2bz) 304 305 ! Debugging section. 306 ABI_CHECK(sum(new%weights) - one < tol6, sjoin("Weights don't sum up to one but to:", ftoa(sum(new%weights)))) 307 308 if (debug /= 0) then 309 do ik_ibz=1,new%nibz 310 if (ik_ibz <= nkibz) then 311 write(std_out,"(a)")sjoin(ktoa(new%ibz(:,ik_ibz)), ktoa(new%ibz(:,ik_ibz) - kibz(:,ik_ibz)), ftoa(new%weights(ik_ibz))) 312 else 313 write(std_out,"(a)")sjoin(ktoa(new%ibz(:,ik_ibz)), "[---]", ftoa(new%weights(ik_ibz))) 314 end if 315 end do 316 end if 317 318 end function lgroup_new
m_sigmaph/lgroup_t [ Types ]
[ Top ] [ m_sigmaph ] [ Types ]
NAME
lgroup_t
FUNCTION
Stores tables associated to the little-group.
SOURCE
65 type, public :: lgroup_t 66 67 integer :: nibz 68 ! Number of points in the IBZ(q) 69 70 integer :: nbz 71 ! Number of points in the full BZ 72 73 integer :: nsym_lg 74 ! Number of operations in the little group. Including "time-reversed" symmetries 75 ! i.e. symmetries S for which -S q = q + G. See also littlegroup_q 76 77 integer :: input_timrev 78 ! 1 if time-reversal symmetry can be used, 0 otherwise. 79 ! NB: This flag refers to the generation of the initial set of k-points. 80 ! One should pay attention when calling other routines in which timrev is required 81 ! Because from Sq = q does not necessarily follow that -Sq = q if q is not on zone-border. 82 ! The operations in G-space stored here already include time-reversal if input_timrev == 1 83 ! so one should call k-point routines with timrev = 0 when using the operations of the little group. 84 85 real(dp) :: point(3) 86 ! The external q-point in reduced coordinates. 87 88 integer,allocatable :: symtab(:,:,:) 89 ! symtab(4, 2, cryst%nsym) 90 ! nsym is the **total** number of spatial symmetries of the system as given by cryst%nsym 91 ! three first numbers define the G vector; 92 ! fourth number is zero if the q-vector is not preserved, is 1 otherwise. 93 ! second index is one without time-reversal symmetry, two with time-reversal symmetry 94 95 integer, allocatable :: symrec_lg(:, :, :) 96 ! symrec_lg(3, 3, nsym_lg) 97 ! Symmetry operations in G-space (including time-reversed operations if any) 98 99 integer, allocatable :: symafm_lg(:) 100 ! symafm_lg(nsym_lg) 101 ! Anti-ferromagnetic character associated to symrec_lg 102 103 integer,allocatable :: bz2ibz_smap(:,:) 104 ! bz2ibz_smap(6, nbz) 105 ! Mapping BZ --> IBZ(q) 106 ! Note that here we used the symmetries of the little group. 107 108 integer, allocatable :: lgsym2glob(:, :) 109 ! lgsym2glob(2, nsym_lg) 110 ! Mapping isym_lg --> [isym, itime] 111 ! where isym is the index of the operation in the global array crystal%symrec 112 ! and itime is 2 if time-reversal T must be included else 1. 113 114 real(dp) :: gmet(3,3) 115 ! Reciprocal space metric in bohr^{-2} 116 117 real(dp),allocatable :: ibz(:,:) 118 ! ibz(3, nibz) 119 ! K-points in the IBZ(q) 120 121 real(dp),allocatable :: weights(:) 122 ! weights(nibz) 123 ! Weights in the IBZ(q), normalized to 1 124 125 contains 126 127 procedure :: findq_ibzk => lgroup_findq_ibzk 128 ! Find the index of the point in the IBZ(k). 129 130 procedure :: find_ibzimage => lgroup_find_ibzimage 131 ! Find the symmetrical image in the IBZ(k) of a qpoint in the BZ. 132 133 procedure :: print => lgroup_print 134 ! Print the object 135 136 procedure :: free => lgroup_free 137 ! Free memory. 138 139 end type lgroup_t