TABLE OF CONTENTS
ABINIT/m_paw_sym [ Modules ]
NAME
m_paw_sym
FUNCTION
This module contains several routines related to the use of symmetries in the PAW approach.
COPYRIGHT
Copyright (C) 2018-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
16 #if defined HAVE_CONFIG_H 17 #include "config.h" 18 #endif 19 20 #include "abi_common.h" 21 22 MODULE m_paw_sym 23 24 use defs_basis 25 use m_abicore 26 use m_errors 27 28 use m_crystal, only : crystal_t 29 use m_pawang, only : pawang_type 30 use m_pawtab, only : pawtab_type 31 use m_pawcprj, only : pawcprj_type, pawcprj_alloc, pawcprj_free, pawcprj_copy 32 use m_bz_mesh, only : kmesh_t 33 34 implicit none 35 36 private 37 38 !public procedures. 39 public :: paw_symcprj ! Symetrize the projections cprj=<n,k|p_i> (p_i=NL PAW projector) - in-place version 40 public :: paw_symcprj_op ! Symetrize the projections cprj=<n,k|p_i> (p_i=NL PAW projector) - out-of-place version 41 42 CONTAINS !========================================================================================
m_paw_sym/paw_symcprj [ Functions ]
[ Top ] [ m_paw_sym ] [ Functions ]
NAME
paw_symcprj
FUNCTION
Symetrize the projections cprj=<n,k|p_i> (p_i=NL PAW projector) - in-place version
INPUTS
ik_ibz=The index of the k-point in the full BZ where the matrix elements have to be symmetrized. nspinor=Number of spinorial components nband_k=Number of bands stored in cprjnk_kibz for this k-point. Cryst<crystal_t>=data type gathering information on unit cell and symmetries. %ntypat=number of type of atoms %natom=number of atoms in the unit cell %typat(natom)=type of each atom %indsym(4,nsym,natom)=indirect indexing array: for each isym,iatom, fourth element is label of atom into which iatom is sent by the INVERSE of the symmetry operation symrel(isym); first three elements are the primitive translations that must be subtracted after the transformation to get back to the original unit cell. Kmesh<kmesh_t>: datatype gathering information on the k-point sampling. %nbz=number of k-points in the full Brillouin zone %nibz=number of k-points in the irreducible wedge %tab(nkbz)=table giving for each k-point in the BZ (array kbz), the corresponding irred. point in the IBZ. i.e k_BZ = (IS) kIBZ where S is one of the symrec operations and I is the inversion or the identity %tabi(nkbz)=for each k-point in the BZ defines whether inversion has to be considered in the relation k_BZ=(IS) k_IBZ (1 => only S; -1 => -S) %tabo(nkbz)= the symmetry operation S that takes k_IBZ to each k_BZ Pawtab(Cryst%ntypat) <type(pawtab_type)>=paw tabulated starting data. Pawang <type(pawang_type)>=paw angular mesh and related data %lmax=Max angular momentum included in the PAW datasets used. mentioned at the second line of the psp file %zarot(2*lmax+1,2*lmax+1,lmax+1,nsym)=coefficients of the transformation of real spherical harmonics under the symmetry operations.
NOTES
Derivatives are not symmetrized.
OUTPUT
SOURCE
87 subroutine paw_symcprj(ik_bz,nspinor,nband_k,Cryst,Kmesh,Pawtab,Pawang,Cprj_bz) 88 89 !Arguments ------------------------------------ 90 !scalars 91 integer,intent(in) :: nspinor,nband_k,ik_bz 92 type(crystal_t),intent(in) :: Cryst 93 type(kmesh_t),intent(in) :: Kmesh 94 type(Pawang_type),intent(in) :: Pawang 95 !arrays 96 type(Pawtab_type),target,intent(in) :: Pawtab(Cryst%ntypat) 97 type(pawcprj_type),intent(inout) :: Cprj_bz(Cryst%natom,nspinor*nband_k) 98 99 !Local variables------------------------------- 100 !scalars 101 integer :: iatom,iat_sym,iband,ibsp_bz 102 integer :: ibsp_ibz,ik_ibz,indexj,ispinor,isym,itim 103 integer :: itypat,jl,jl0,jlmn,jln,jln0,jlpm,jm,jn,lmax,mm,ncpgr 104 real(dp) :: arg,wtk 105 logical :: isirred 106 !arrays 107 integer :: r0(3),nlmn_atom(Cryst%natom) 108 integer,ABI_CONTIGUOUS pointer :: indlmn(:,:) 109 real(dp) :: dum(2,nspinor),kbz(3),kirr(3),phase(2),swp(2),tmp(2,nspinor) 110 real(dp),allocatable :: DS_mmpl(:,:,:) 111 type(pawcprj_type) :: Cprjnk_kibz(Cryst%natom,nspinor*nband_k) 112 113 ! ********************************************************************* 114 115 ncpgr = Cprj_bz(1,1)%ncpgr 116 ABI_CHECK(ncpgr==0,"Derivatives of cprj are not coded") 117 118 !Get the index of the IBZ image associated to the BZ k-point ik_bz and related simmetry. 119 call Kmesh%get_BZ_item(ik_bz,kbz,ik_ibz,isym,itim,isirred=isirred) 120 121 if (isirred) RETURN ! It is a point in the IBZ, Symmetrization is not needed. 122 ! 123 !The corresponding point kirr in the IBZ. 124 call kmesh%get_IBZ_item(ik_ibz,kirr,wtk) 125 126 !Local copy. 127 do iatom=1,Cryst%natom 128 nlmn_atom(iatom)=Pawtab(Cryst%typat(iatom))%lmn_size 129 end do 130 131 call pawcprj_alloc(Cprjnk_kibz,ncpgr,nlmn_atom) 132 call pawcprj_copy(Cprj_bz,Cprjnk_kibz) 133 ! 134 !=== DS_mmpl is the rotation matrix for real spherical harmonics associated to symrec(:,:,isym) === 135 !* Note the convention used by Blanco in Eq. 27 : DS_mmp multiply spherical harmonics as row vectors 136 lmax=Pawang%l_max-1 ! l_max is Max l+1 137 ABI_MALLOC(DS_mmpl,(2*lmax+1,2*lmax+1,lmax+1)) 138 DS_mmpl=Pawang%zarot(:,:,:,isym) 139 ! 140 !=========================================== 141 !==== Loop over atoms to be symmetrized ==== 142 !=========================================== 143 do iatom=1,Cryst%natom 144 itypat=Cryst%typat(iatom) 145 iat_sym=Cryst%indsym(4,isym,iatom) 146 indlmn => Pawtab(itypat)%indlmn 147 r0=Cryst%indsym(1:3,isym,iatom) ! R^{-1} (xred(:,iatom)-tnons) = xred(:,iat_sym) + r0. 148 arg=two_pi*dot_product(kirr,r0) 149 phase(1)=COS(arg) 150 phase(2)=SIN(arg) 151 ! 152 ! Loop over the (jl,jm,jn) components to be symmetrized. 153 jl0=-1; jln0=-1; indexj=1 154 do jlmn=1,Pawtab(itypat)%lmn_size 155 jl =indlmn(1,jlmn) 156 jm =indlmn(2,jlmn) 157 jn =indlmn(3,jlmn) 158 jln =indlmn(5,jlmn) 159 jlpm=1+jl+jm 160 if (jln/=jln0) indexj=indexj+2*jl0+1 161 ! 162 ! === For each band, calculate contribution due to rotated real spherical harmonics === 163 ! FIXME check this expression; according to Blanco I should have D(S^-1} but it seems D(S) is correct 164 ! Recheck spinorial case, presently is wrong 165 ibsp_ibz=0 166 ibsp_bz=0 167 do iband=1,nband_k 168 169 tmp(:,:)=zero 170 do ispinor=1,nspinor 171 ibsp_ibz=ibsp_ibz+1 172 do mm=1,2*jl+1 173 tmp(1,ispinor)=tmp(1,ispinor)+DS_mmpl(mm,jlpm,jl+1)*Cprjnk_kibz(iat_sym,ibsp_ibz)%cp(1,indexj+mm) 174 tmp(2,ispinor)=tmp(2,ispinor)+DS_mmpl(mm,jlpm,jl+1)*Cprjnk_kibz(iat_sym,ibsp_ibz)%cp(2,indexj+mm) 175 end do 176 end do !ispinor 177 ! 178 ! * Apply the phase to account if the symmetric atom belongs to a different unit cell. 179 do ispinor=1,nspinor 180 dum(1,ispinor)=tmp(1,ispinor)*phase(1)-tmp(2,ispinor)*phase(2) 181 dum(2,ispinor)=tmp(1,ispinor)*phase(2)+tmp(2,ispinor)*phase(1) 182 end do 183 ! 184 ! * If required, apply time-reversal symmetry to retrieve the correct point in the BZ. 185 if (itim==2) then 186 if (nspinor==1) then 187 dum(2,1)=-dum(2,1) 188 else if (nspinor==2) then ! TODO rotate wavefunction in spinor space. 189 swp(:)=dum(:,1) 190 dum(1,1)= dum(1,2) 191 dum(2,1)=-dum(2,2) 192 dum(1,2)=-swp(1) 193 dum(2,2)= swp(2) 194 end if 195 end if 196 ! 197 ! ==== Save values ==== 198 do ispinor=1,nspinor 199 ibsp_bz=ibsp_bz+1 200 Cprj_bz(iatom,ibsp_bz)%cp(1,jlmn)=dum(1,ispinor) 201 Cprj_bz(iatom,ibsp_bz)%cp(2,jlmn)=dum(2,ispinor) 202 end do 203 end do !iband 204 205 jl0=jl; jln0=jln 206 end do !jlmn 207 end do !iatom 208 209 call pawcprj_free(Cprjnk_kibz) 210 ABI_FREE(DS_mmpl) 211 212 end subroutine paw_symcprj
m_paw_sym/paw_symcprj_op [ Functions ]
[ Top ] [ m_paw_sym ] [ Functions ]
NAME
paw_symcprj_op
FUNCTION
Symetrize the projections cprj=<n,k|p_i> (p_i=NL PAW projector) - in-place version
INPUTS
ik_ibz=The index of the k-point in the full BZ where the matrix elements have to be symmetrized. nspinor=Number of spinorial components nband_k=Number of bands stored in cprjnk_kibz for this k-point. Cryst<crystal_t>=data type gathering information on unit cell and symmetries. %ntypat=number of type of atoms %natom=number of atoms in the unit cell %typat(natom)=type of each atom %indsym(4,nsym,natom)=indirect indexing array: for each isym,iatom, fourth element is label of atom into which iatom is sent by the INVERSE of the symmetry operation symrel(isym); first three elements are the primitive translations that must be subtracted after the transformation to get back to the original unit cell. Kmesh<kmesh_t>: datatype gathering information on the k-point sampling. %nbz=number of k-points in the full Brillouin zone %nibz=number of k-points in the irreducible wedge %tab(nkbz)=table giving for each k-point in the BZ (array kbz), the corresponding irred. point in the IBZ. i.e k_BZ = (IS) kIBZ where S is one of the symrec operations and I is the inversion or the identity %tabi(nkbz)=for each k-point in the BZ defines whether inversion has to be considered in the relation k_BZ=(IS) k_IBZ (1 => only S; -1 => -S) %tabo(nkbz)= the symmetry operation S that takes k_IBZ to each k_BZ Pawtab(Cryst%ntypat) <type(pawtab_type)>=paw tabulated starting data. Pawang <type(pawang_type)>=paw angular mesh and related data %lmax=Max angular momentum included in the PAW datasets used. mentioned at the second line of the psp file %zarot(2*lmax+1,2*lmax+1,lmax+1,nsym)=coefficients of the transformation of real spherical harmonics under the symmetry operations. in_Cprj(Cryst%natom,nspinor*nband_k)<pawcprj_type>=Input cprj
OUTPUT
out_Cprj(Cryst%natom,nspinor*nband_k)<pawcprj_type>=Symmetrized cprj matrix elements.
NOTES
Derivatives are not symmetrized.
SOURCE
259 subroutine paw_symcprj_op(ik_bz,nspinor,nband_k,Cryst,Kmesh,Pawtab,Pawang,in_Cprj,out_Cprj) 260 261 !Arguments ------------------------------------ 262 !scalars 263 integer,intent(in) :: nspinor,nband_k,ik_bz 264 type(crystal_t),intent(in) :: Cryst 265 type(kmesh_t),intent(in) :: Kmesh 266 type(Pawang_type),intent(in) :: Pawang 267 !arrays 268 type(Pawtab_type),target,intent(in) :: Pawtab(Cryst%ntypat) 269 type(pawcprj_type),intent(in) :: in_Cprj(Cryst%natom,nspinor*nband_k) 270 type(pawcprj_type),intent(inout) :: out_Cprj(Cryst%natom,nspinor*nband_k) !vz_i 271 272 !Local variables------------------------------- 273 !scalars 274 integer :: iatom,iat_sym,iband,ibsp_bz 275 integer :: ibsp_ibz,ik_ibz,indexj,ispinor,isym,itim 276 integer :: itypat,jl,jl0,jlmn,jln,jln0,jlpm,jm,jn,lmax,mm,ncpgr 277 real(dp) :: arg,wtk 278 logical :: isirred 279 !arrays 280 integer :: r0(3) !,nlmn_atom(Cryst%natom) 281 integer,ABI_CONTIGUOUS pointer :: indlmn(:,:) 282 real(dp) :: dum(2,nspinor),kbz(3),kirr(3),phase(2),swp(2),tmp(2,nspinor) 283 real(dp),allocatable :: DS_mmpl(:,:,:) 284 285 ! ********************************************************************* 286 287 ncpgr = in_Cprj(1,1)%ncpgr 288 ABI_CHECK(ncpgr==0,"Derivatives of cprj are not coded") 289 290 !Get the index of the IBZ image associated to the BZ k-point ik_bz and related simmetry. 291 call Kmesh%get_BZ_item(ik_bz,kbz,ik_ibz,isym,itim,isirred=isirred) 292 293 if (isirred) then ! It is a point in the IBZ, Symmetrization is not needed. 294 call pawcprj_copy(in_Cprj,out_Cprj) 295 RETURN 296 end if 297 ! 298 !The corresponding point kirr in the IBZ. 299 call Kmesh%get_IBZ_item(ik_ibz,kirr,wtk) 300 ! 301 !=== DS_mmpl is the rotation matrix for real spherical harmonics associated to symrec(:,:,isym) === 302 !* Note the convention used by Blanco in Eq. 27 : DS_mmp multiply spherical harmonics as row vectors 303 lmax=Pawang%l_max-1 ! l_max is Max l+1 304 ABI_MALLOC(DS_mmpl,(2*lmax+1,2*lmax+1,lmax+1)) 305 DS_mmpl=Pawang%zarot(:,:,:,isym) 306 307 !Local copy. 308 !do iatom=1,Cryst%natom 309 !nlmn_atom(iatom)=Pawtab(Cryst%typat(iatom))%lmn_size 310 !end do 311 !call pawcprj_alloc(out_Cprj,ncpgr,nlmn_atom) 312 ! 313 !=========================================== 314 !==== Loop over atoms to be symmetrized ==== 315 !=========================================== 316 do iatom=1,Cryst%natom 317 itypat=Cryst%typat(iatom) 318 iat_sym=Cryst%indsym(4,isym,iatom) 319 indlmn => Pawtab(itypat)%indlmn 320 r0=Cryst%indsym(1:3,isym,iatom) ! R^{-1} (xred(:,iatom)-tnons) = xred(:,iat_sym) + r0. 321 arg=two_pi*dot_product(kirr,r0) 322 phase(1)=COS(arg) 323 phase(2)=SIN(arg) 324 ! 325 ! Loop over the (jl,jm,jn) components to be symmetrized. 326 jl0=-1; jln0=-1; indexj=1 327 do jlmn=1,Pawtab(itypat)%lmn_size 328 jl =indlmn(1,jlmn) 329 jm =indlmn(2,jlmn) 330 jn =indlmn(3,jlmn) 331 jln =indlmn(5,jlmn) 332 jlpm=1+jl+jm 333 if (jln/=jln0) indexj=indexj+2*jl0+1 334 ! 335 ! === For each band, calculate contribution due to rotated real spherical harmonics === 336 ! FIXME check this expression; according to Blanco I should have D(S^-1} but it seems D(S) is correct 337 ! Recheck spinorial case, presently is wrong 338 ibsp_ibz=0 339 ibsp_bz=0 340 do iband=1,nband_k 341 342 tmp(:,:)=zero 343 do ispinor=1,nspinor 344 ibsp_ibz=ibsp_ibz+1 345 do mm=1,2*jl+1 346 tmp(1,ispinor)=tmp(1,ispinor)+DS_mmpl(mm,jlpm,jl+1)*in_Cprj(iat_sym,ibsp_ibz)%cp(1,indexj+mm) 347 tmp(2,ispinor)=tmp(2,ispinor)+DS_mmpl(mm,jlpm,jl+1)*in_Cprj(iat_sym,ibsp_ibz)%cp(2,indexj+mm) 348 end do 349 end do !ispinor 350 ! 351 ! * Apply the phase to account if the symmetric atom belongs to a different unit cell. 352 do ispinor=1,nspinor 353 dum(1,ispinor)=tmp(1,ispinor)*phase(1)-tmp(2,ispinor)*phase(2) 354 dum(2,ispinor)=tmp(1,ispinor)*phase(2)+tmp(2,ispinor)*phase(1) 355 end do 356 ! 357 ! * If required, apply time-reversal symmetry to retrieve the correct point in the BZ. 358 if (itim==2) then 359 if (nspinor==1) then 360 dum(2,1)=-dum(2,1) 361 else if (nspinor==2) then ! TODO rotate wavefunction in spinor space. 362 swp(:)=dum(:,1) 363 dum(1,1)= dum(1,2) 364 dum(2,1)=-dum(2,2) 365 dum(1,2)=-swp(1) 366 dum(2,2)= swp(2) 367 end if 368 end if 369 ! 370 ! ==== Save values ==== 371 do ispinor=1,nspinor 372 ibsp_bz=ibsp_bz+1 373 out_Cprj(iatom,ibsp_bz)%cp(1,jlmn)=dum(1,ispinor) 374 out_Cprj(iatom,ibsp_bz)%cp(2,jlmn)=dum(2,ispinor) 375 end do 376 end do !iband 377 378 jl0=jl; jln0=jln 379 end do !jlmn 380 end do !iatom 381 382 ABI_FREE(DS_mmpl) 383 384 end subroutine paw_symcprj_op