TABLE OF CONTENTS
ABINIT/paw_symcprj [ Functions ]
NAME
paw_symcprj
FUNCTION
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
PARENTS
calc_sigc_me,calc_sigx_me,cchi0,cchi0q0,cchi0q0_intraband,cohsex_me debug_tools,m_shirley,prep_calc_ucrpa
CHILDREN
get_bz_item,get_ibz_item,pawcprj_copy
SOURCE
48 #if defined HAVE_CONFIG_H 49 #include "config.h" 50 #endif 51 52 #include "abi_common.h" 53 54 55 subroutine paw_symcprj(ik_bz,nspinor,nband_k,Cryst,Kmesh,Pawtab,Pawang,Cprj_bz) 56 57 use defs_basis 58 use m_profiling_abi 59 use m_errors 60 61 use m_crystal, only : crystal_t 62 use m_pawang, only : pawang_type 63 use m_pawtab, only : pawtab_type 64 use m_pawcprj, only : pawcprj_type, pawcprj_alloc, pawcprj_free, pawcprj_copy 65 use m_bz_mesh, only : kmesh_t, get_ibz_item, get_bz_item 66 67 !This section has been created automatically by the script Abilint (TD). 68 !Do not modify the following lines by hand. 69 #undef ABI_FUNC 70 #define ABI_FUNC 'paw_symcprj' 71 !End of the abilint section 72 73 implicit none 74 75 !Arguments ------------------------------------ 76 !scalars 77 integer,intent(in) :: nspinor,nband_k,ik_bz 78 type(crystal_t),intent(in) :: Cryst 79 type(kmesh_t),intent(in) :: Kmesh 80 type(Pawang_type),intent(in) :: Pawang 81 !arrays 82 type(Pawtab_type),target,intent(in) :: Pawtab(Cryst%ntypat) 83 type(pawcprj_type),intent(inout) :: Cprj_bz(Cryst%natom,nspinor*nband_k) 84 85 !Local variables------------------------------- 86 !scalars 87 integer :: iatom,iat_sym,iband,ibsp_bz 88 integer :: ibsp_ibz,ik_ibz,indexj,ispinor,isym,itim 89 integer :: itypat,jl,jl0,jlmn,jln,jln0,jlpm,jm,jn,lmax,mm,ncpgr 90 real(dp) :: arg,wtk 91 logical :: isirred 92 !arrays 93 integer :: r0(3),nlmn_atom(Cryst%natom) 94 integer,ABI_CONTIGUOUS pointer :: indlmn(:,:) 95 real(dp) :: dum(2,nspinor),kbz(3),kirr(3),phase(2),swp(2),tmp(2,nspinor) 96 real(dp),allocatable :: DS_mmpl(:,:,:) 97 type(pawcprj_type) :: Cprjnk_kibz(Cryst%natom,nspinor*nband_k) 98 99 ! ********************************************************************* 100 101 ncpgr = Cprj_bz(1,1)%ncpgr 102 ABI_CHECK(ncpgr==0,"Derivatived of cprj are not coded") 103 104 !Get the index of the IBZ image associated to the BZ k-point ik_bz and related simmetry. 105 call get_BZ_item(Kmesh,ik_bz,kbz,ik_ibz,isym,itim,isirred=isirred) 106 107 if (isirred) RETURN ! It is a point in the IBZ, Symmetrization is not needed. 108 ! 109 !The corresponding point kirr in the IBZ. 110 call get_IBZ_item(Kmesh,ik_ibz,kirr,wtk) 111 112 !Local copy. 113 do iatom=1,Cryst%natom 114 nlmn_atom(iatom)=Pawtab(Cryst%typat(iatom))%lmn_size 115 end do 116 117 call pawcprj_alloc(Cprjnk_kibz,ncpgr,nlmn_atom) 118 call pawcprj_copy(Cprj_bz,Cprjnk_kibz) 119 ! 120 !=== DS_mmpl is the rotation matrix for real spherical harmonics associated to symrec(:,:,isym) === 121 !* Note the convention used by Blanco in Eq. 27 : DS_mmp multiply spherical harmonics as row vectors 122 lmax=Pawang%l_max-1 ! l_max is Max l+1 123 ABI_ALLOCATE(DS_mmpl,(2*lmax+1,2*lmax+1,lmax+1)) 124 DS_mmpl=Pawang%zarot(:,:,:,isym) 125 ! 126 !=========================================== 127 !==== Loop over atoms to be symmetrized ==== 128 !=========================================== 129 do iatom=1,Cryst%natom 130 itypat=Cryst%typat(iatom) 131 iat_sym=Cryst%indsym(4,isym,iatom) 132 indlmn => Pawtab(itypat)%indlmn 133 r0=Cryst%indsym(1:3,isym,iatom) ! R^{-1} (xred(:,iatom)-tnons) = xred(:,iat_sym) + r0. 134 arg=two_pi*dot_product(kirr,r0) 135 phase(1)=COS(arg) 136 phase(2)=SIN(arg) 137 ! 138 ! Loop over the (jl,jm,jn) components to be symmetrized. 139 jl0=-1; jln0=-1; indexj=1 140 do jlmn=1,Pawtab(itypat)%lmn_size 141 jl =indlmn(1,jlmn) 142 jm =indlmn(2,jlmn) 143 jn =indlmn(3,jlmn) 144 jln =indlmn(5,jlmn) 145 jlpm=1+jl+jm 146 if (jln/=jln0) indexj=indexj+2*jl0+1 147 ! 148 ! === For each band, calculate contribution due to rotated real spherical harmonics === 149 ! FIXME check this expression; according to Blanco I should have D(S^-1} but it seems D(S) is correct 150 ! Recheck spinorial case, presently is wrong 151 ibsp_ibz=0 152 ibsp_bz=0 153 do iband=1,nband_k 154 155 tmp(:,:)=zero 156 do ispinor=1,nspinor 157 ibsp_ibz=ibsp_ibz+1 158 do mm=1,2*jl+1 159 tmp(1,ispinor)=tmp(1,ispinor)+DS_mmpl(mm,jlpm,jl+1)*Cprjnk_kibz(iat_sym,ibsp_ibz)%cp(1,indexj+mm) 160 tmp(2,ispinor)=tmp(2,ispinor)+DS_mmpl(mm,jlpm,jl+1)*Cprjnk_kibz(iat_sym,ibsp_ibz)%cp(2,indexj+mm) 161 end do 162 end do !ispinor 163 ! 164 ! * Apply the phase to account if the symmetric atom belongs to a different unit cell. 165 do ispinor=1,nspinor 166 dum(1,ispinor)=tmp(1,ispinor)*phase(1)-tmp(2,ispinor)*phase(2) 167 dum(2,ispinor)=tmp(1,ispinor)*phase(2)+tmp(2,ispinor)*phase(1) 168 end do 169 ! 170 ! * If required, apply time-reversal symmetry to retrieve the correct point in the BZ. 171 if (itim==2) then 172 if (nspinor==1) then 173 dum(2,1)=-dum(2,1) 174 else if (nspinor==2) then ! TODO rotate wavefunction in spinor space. 175 swp(:)=dum(:,1) 176 dum(1,1)= dum(1,2) 177 dum(2,1)=-dum(2,2) 178 dum(1,2)=-swp(1) 179 dum(2,2)= swp(2) 180 end if 181 end if 182 ! 183 ! ==== Save values ==== 184 do ispinor=1,nspinor 185 ibsp_bz=ibsp_bz+1 186 Cprj_bz(iatom,ibsp_bz)%cp(1,jlmn)=dum(1,ispinor) 187 Cprj_bz(iatom,ibsp_bz)%cp(2,jlmn)=dum(2,ispinor) 188 end do 189 end do !iband 190 191 jl0=jl; jln0=jln 192 end do !jlmn 193 end do !iatom 194 195 call pawcprj_free(Cprjnk_kibz) 196 ABI_DEALLOCATE(DS_mmpl) 197 198 end subroutine paw_symcprj
ABINIT/paw_symcprj_op [ Functions ]
NAME
paw_symcprj_op
FUNCTION
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.
PARENTS
exc_build_block
CHILDREN
get_bz_item,get_ibz_item,pawcprj_copy
SOURCE
250 subroutine paw_symcprj_op(ik_bz,nspinor,nband_k,Cryst,Kmesh,Pawtab,Pawang,in_Cprj,out_Cprj) 251 252 use defs_basis 253 use m_profiling_abi 254 use m_errors 255 256 use m_crystal, only : crystal_t 257 use m_bz_mesh, only : kmesh_t, get_ibz_item, get_bz_item 258 use m_pawang, only : pawang_type 259 use m_pawtab, only : pawtab_type 260 use m_pawcprj, only : pawcprj_type, pawcprj_alloc, pawcprj_free, pawcprj_copy 261 262 !This section has been created automatically by the script Abilint (TD). 263 !Do not modify the following lines by hand. 264 #undef ABI_FUNC 265 #define ABI_FUNC 'paw_symcprj_op' 266 !End of the abilint section 267 268 implicit none 269 270 !Arguments ------------------------------------ 271 !scalars 272 integer,intent(in) :: nspinor,nband_k,ik_bz 273 type(crystal_t),intent(in) :: Cryst 274 type(kmesh_t),intent(in) :: Kmesh 275 type(Pawang_type),intent(in) :: Pawang 276 !arrays 277 type(Pawtab_type),target,intent(in) :: Pawtab(Cryst%ntypat) 278 type(pawcprj_type),intent(in) :: in_Cprj(Cryst%natom,nspinor*nband_k) 279 type(pawcprj_type),intent(inout) :: out_Cprj(Cryst%natom,nspinor*nband_k) !vz_i 280 281 !Local variables------------------------------- 282 !scalars 283 integer :: iatom,iat_sym,iband,ibsp_bz 284 integer :: ibsp_ibz,ik_ibz,indexj,ispinor,isym,itim 285 integer :: itypat,jl,jl0,jlmn,jln,jln0,jlpm,jm,jn,lmax,mm,ncpgr 286 real(dp) :: arg,wtk 287 logical :: isirred 288 !arrays 289 integer :: r0(3) !,nlmn_atom(Cryst%natom) 290 integer,ABI_CONTIGUOUS pointer :: indlmn(:,:) 291 real(dp) :: dum(2,nspinor),kbz(3),kirr(3),phase(2),swp(2),tmp(2,nspinor) 292 real(dp),allocatable :: DS_mmpl(:,:,:) 293 294 ! ********************************************************************* 295 296 ncpgr = in_Cprj(1,1)%ncpgr 297 ABI_CHECK(ncpgr==0,"Derivatived of cprj are not coded") 298 299 !Get the index of the IBZ image associated to the BZ k-point ik_bz and related simmetry. 300 call get_BZ_item(Kmesh,ik_bz,kbz,ik_ibz,isym,itim,isirred=isirred) 301 302 if (isirred) then ! It is a point in the IBZ, Symmetrization is not needed. 303 call pawcprj_copy(in_Cprj,out_Cprj) 304 RETURN 305 end if 306 ! 307 !The corresponding point kirr in the IBZ. 308 call get_IBZ_item(Kmesh,ik_ibz,kirr,wtk) 309 ! 310 !=== DS_mmpl is the rotation matrix for real spherical harmonics associated to symrec(:,:,isym) === 311 !* Note the convention used by Blanco in Eq. 27 : DS_mmp multiply spherical harmonics as row vectors 312 lmax=Pawang%l_max-1 ! l_max is Max l+1 313 ABI_ALLOCATE(DS_mmpl,(2*lmax+1,2*lmax+1,lmax+1)) 314 DS_mmpl=Pawang%zarot(:,:,:,isym) 315 316 !Local copy. 317 !do iatom=1,Cryst%natom 318 !nlmn_atom(iatom)=Pawtab(Cryst%typat(iatom))%lmn_size 319 !end do 320 !call pawcprj_alloc(out_Cprj,ncpgr,nlmn_atom) 321 ! 322 !=========================================== 323 !==== Loop over atoms to be symmetrized ==== 324 !=========================================== 325 do iatom=1,Cryst%natom 326 itypat=Cryst%typat(iatom) 327 iat_sym=Cryst%indsym(4,isym,iatom) 328 indlmn => Pawtab(itypat)%indlmn 329 r0=Cryst%indsym(1:3,isym,iatom) ! R^{-1} (xred(:,iatom)-tnons) = xred(:,iat_sym) + r0. 330 arg=two_pi*dot_product(kirr,r0) 331 phase(1)=COS(arg) 332 phase(2)=SIN(arg) 333 ! 334 ! Loop over the (jl,jm,jn) components to be symmetrized. 335 jl0=-1; jln0=-1; indexj=1 336 do jlmn=1,Pawtab(itypat)%lmn_size 337 jl =indlmn(1,jlmn) 338 jm =indlmn(2,jlmn) 339 jn =indlmn(3,jlmn) 340 jln =indlmn(5,jlmn) 341 jlpm=1+jl+jm 342 if (jln/=jln0) indexj=indexj+2*jl0+1 343 ! 344 ! === For each band, calculate contribution due to rotated real spherical harmonics === 345 ! FIXME check this expression; according to Blanco I should have D(S^-1} but it seems D(S) is correct 346 ! Recheck spinorial case, presently is wrong 347 ibsp_ibz=0 348 ibsp_bz=0 349 do iband=1,nband_k 350 351 tmp(:,:)=zero 352 do ispinor=1,nspinor 353 ibsp_ibz=ibsp_ibz+1 354 do mm=1,2*jl+1 355 tmp(1,ispinor)=tmp(1,ispinor)+DS_mmpl(mm,jlpm,jl+1)*in_Cprj(iat_sym,ibsp_ibz)%cp(1,indexj+mm) 356 tmp(2,ispinor)=tmp(2,ispinor)+DS_mmpl(mm,jlpm,jl+1)*in_Cprj(iat_sym,ibsp_ibz)%cp(2,indexj+mm) 357 end do 358 end do !ispinor 359 ! 360 ! * Apply the phase to account if the symmetric atom belongs to a different unit cell. 361 do ispinor=1,nspinor 362 dum(1,ispinor)=tmp(1,ispinor)*phase(1)-tmp(2,ispinor)*phase(2) 363 dum(2,ispinor)=tmp(1,ispinor)*phase(2)+tmp(2,ispinor)*phase(1) 364 end do 365 ! 366 ! * If required, apply time-reversal symmetry to retrieve the correct point in the BZ. 367 if (itim==2) then 368 if (nspinor==1) then 369 dum(2,1)=-dum(2,1) 370 else if (nspinor==2) then ! TODO rotate wavefunction in spinor space. 371 swp(:)=dum(:,1) 372 dum(1,1)= dum(1,2) 373 dum(2,1)=-dum(2,2) 374 dum(1,2)=-swp(1) 375 dum(2,2)= swp(2) 376 end if 377 end if 378 ! 379 ! ==== Save values ==== 380 do ispinor=1,nspinor 381 ibsp_bz=ibsp_bz+1 382 out_Cprj(iatom,ibsp_bz)%cp(1,jlmn)=dum(1,ispinor) 383 out_Cprj(iatom,ibsp_bz)%cp(2,jlmn)=dum(2,ispinor) 384 end do 385 end do !iband 386 387 jl0=jl; jln0=jln 388 end do !jlmn 389 end do !iatom 390 391 ABI_DEALLOCATE(DS_mmpl) 392 393 end subroutine paw_symcprj_op