TABLE OF CONTENTS


ABINIT/m_paw_sym [ Modules ]

[ Top ] [ 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-2018 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

17 #if defined HAVE_CONFIG_H
18 #include "config.h"
19 #endif
20 
21 #include "abi_common.h"
22 
23 MODULE m_paw_sym
24 
25  use defs_basis
26  use m_abicore
27  use m_errors
28 
29  use m_crystal,   only : crystal_t
30  use m_pawang,    only : pawang_type
31  use m_pawtab,    only : pawtab_type
32  use m_pawcprj,   only : pawcprj_type,pawcprj_alloc,pawcprj_free,pawcprj_copy
33  use m_bz_mesh,   only : kmesh_t,get_ibz_item,get_bz_item
34 
35  implicit none
36 
37  private
38 
39 !public procedures.
40  public :: paw_symcprj    ! Symetrize the projections cprj=<n,k|p_i> (p_i=NL PAW projector) - in-place version
41  public :: paw_symcprj_op ! Symetrize the projections cprj=<n,k|p_i> (p_i=NL PAW projector) - out-of-place version
42 
43 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

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

 95 subroutine paw_symcprj(ik_bz,nspinor,nband_k,Cryst,Kmesh,Pawtab,Pawang,Cprj_bz)
 96 
 97 
 98 !This section has been created automatically by the script Abilint (TD).
 99 !Do not modify the following lines by hand.
100 #undef ABI_FUNC
101 #define ABI_FUNC 'paw_symcprj'
102 !End of the abilint section
103 
104  implicit none
105 
106 !Arguments ------------------------------------
107 !scalars
108  integer,intent(in) :: nspinor,nband_k,ik_bz
109  type(crystal_t),intent(in) :: Cryst
110  type(kmesh_t),intent(in) :: Kmesh
111  type(Pawang_type),intent(in) :: Pawang
112 !arrays
113  type(Pawtab_type),target,intent(in) :: Pawtab(Cryst%ntypat)
114  type(pawcprj_type),intent(inout) :: Cprj_bz(Cryst%natom,nspinor*nband_k)
115 
116 !Local variables-------------------------------
117 !scalars
118  integer :: iatom,iat_sym,iband,ibsp_bz
119  integer :: ibsp_ibz,ik_ibz,indexj,ispinor,isym,itim
120  integer :: itypat,jl,jl0,jlmn,jln,jln0,jlpm,jm,jn,lmax,mm,ncpgr
121  real(dp) :: arg,wtk
122  logical :: isirred
123 !arrays
124  integer :: r0(3),nlmn_atom(Cryst%natom)
125  integer,ABI_CONTIGUOUS pointer :: indlmn(:,:)
126  real(dp) :: dum(2,nspinor),kbz(3),kirr(3),phase(2),swp(2),tmp(2,nspinor)
127  real(dp),allocatable :: DS_mmpl(:,:,:)
128  type(pawcprj_type) :: Cprjnk_kibz(Cryst%natom,nspinor*nband_k)
129 
130 ! *********************************************************************
131 
132  ncpgr = Cprj_bz(1,1)%ncpgr
133  ABI_CHECK(ncpgr==0,"Derivatives of cprj are not coded")
134 
135 !Get the index of the IBZ image associated to the BZ k-point ik_bz and related simmetry.
136  call get_BZ_item(Kmesh,ik_bz,kbz,ik_ibz,isym,itim,isirred=isirred)
137 
138  if (isirred) RETURN  ! It is a point in the IBZ, Symmetrization is not needed.
139 !
140 !The corresponding point kirr in the IBZ.
141  call get_IBZ_item(Kmesh,ik_ibz,kirr,wtk)
142 
143 !Local copy.
144  do iatom=1,Cryst%natom
145    nlmn_atom(iatom)=Pawtab(Cryst%typat(iatom))%lmn_size
146  end do
147  
148  call pawcprj_alloc(Cprjnk_kibz,ncpgr,nlmn_atom)
149  call pawcprj_copy(Cprj_bz,Cprjnk_kibz)
150 !
151 !=== DS_mmpl is the rotation matrix for real spherical harmonics associated to symrec(:,:,isym) ===
152 !* Note the convention used by Blanco in Eq. 27 : DS_mmp multiply spherical harmonics as row vectors
153  lmax=Pawang%l_max-1 ! l_max is Max l+1
154  ABI_ALLOCATE(DS_mmpl,(2*lmax+1,2*lmax+1,lmax+1))
155  DS_mmpl=Pawang%zarot(:,:,:,isym)
156 !
157 !===========================================
158 !==== Loop over atoms to be symmetrized ====
159 !===========================================
160  do iatom=1,Cryst%natom
161    itypat=Cryst%typat(iatom)
162    iat_sym=Cryst%indsym(4,isym,iatom)
163    indlmn => Pawtab(itypat)%indlmn
164    r0=Cryst%indsym(1:3,isym,iatom) ! R^{-1} (xred(:,iatom)-tnons) = xred(:,iat_sym) + r0.
165    arg=two_pi*dot_product(kirr,r0)
166    phase(1)=COS(arg)
167    phase(2)=SIN(arg)
168 !
169 !  Loop over the (jl,jm,jn) components to be symmetrized.
170    jl0=-1; jln0=-1; indexj=1
171    do jlmn=1,Pawtab(itypat)%lmn_size
172      jl  =indlmn(1,jlmn)
173      jm  =indlmn(2,jlmn)
174      jn  =indlmn(3,jlmn)
175      jln =indlmn(5,jlmn)
176      jlpm=1+jl+jm
177      if (jln/=jln0) indexj=indexj+2*jl0+1
178 !
179 !    === For each band, calculate contribution due to rotated real spherical harmonics ===
180 !    FIXME check this expression; according to Blanco I should have D(S^-1} but it seems D(S) is correct
181 !    Recheck spinorial case, presently is wrong
182      ibsp_ibz=0
183      ibsp_bz=0
184      do iband=1,nband_k
185 
186        tmp(:,:)=zero
187        do ispinor=1,nspinor
188          ibsp_ibz=ibsp_ibz+1
189          do mm=1,2*jl+1
190            tmp(1,ispinor)=tmp(1,ispinor)+DS_mmpl(mm,jlpm,jl+1)*Cprjnk_kibz(iat_sym,ibsp_ibz)%cp(1,indexj+mm)
191            tmp(2,ispinor)=tmp(2,ispinor)+DS_mmpl(mm,jlpm,jl+1)*Cprjnk_kibz(iat_sym,ibsp_ibz)%cp(2,indexj+mm)
192          end do
193        end do !ispinor
194 !
195 !      * Apply the phase to account if the symmetric atom belongs to a different unit cell.
196        do ispinor=1,nspinor
197          dum(1,ispinor)=tmp(1,ispinor)*phase(1)-tmp(2,ispinor)*phase(2)
198          dum(2,ispinor)=tmp(1,ispinor)*phase(2)+tmp(2,ispinor)*phase(1)
199        end do
200 !
201 !      * If required, apply time-reversal symmetry to retrieve the correct point in the BZ.
202        if (itim==2) then
203          if (nspinor==1) then
204            dum(2,1)=-dum(2,1)
205          else if (nspinor==2) then ! TODO rotate wavefunction in spinor space.
206            swp(:)=dum(:,1)
207            dum(1,1)= dum(1,2)
208            dum(2,1)=-dum(2,2)
209            dum(1,2)=-swp(1)
210            dum(2,2)= swp(2)
211          end if
212        end if
213 !
214 !      ==== Save values ====
215        do ispinor=1,nspinor
216          ibsp_bz=ibsp_bz+1
217          Cprj_bz(iatom,ibsp_bz)%cp(1,jlmn)=dum(1,ispinor)
218          Cprj_bz(iatom,ibsp_bz)%cp(2,jlmn)=dum(2,ispinor)
219        end do
220      end do !iband
221 
222      jl0=jl; jln0=jln
223    end do !jlmn
224  end do !iatom
225 
226  call pawcprj_free(Cprjnk_kibz)
227  ABI_DEALLOCATE(DS_mmpl)
228 
229 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.

PARENTS

      exc_build_block

CHILDREN

      get_bz_item,get_ibz_item,pawcprj_copy

SOURCE

282 subroutine paw_symcprj_op(ik_bz,nspinor,nband_k,Cryst,Kmesh,Pawtab,Pawang,in_Cprj,out_Cprj)
283 
284 
285 !This section has been created automatically by the script Abilint (TD).
286 !Do not modify the following lines by hand.
287 #undef ABI_FUNC
288 #define ABI_FUNC 'paw_symcprj_op'
289 !End of the abilint section
290 
291  implicit none
292 
293 !Arguments ------------------------------------
294 !scalars
295  integer,intent(in) :: nspinor,nband_k,ik_bz
296  type(crystal_t),intent(in) :: Cryst
297  type(kmesh_t),intent(in) :: Kmesh
298  type(Pawang_type),intent(in) :: Pawang
299 !arrays
300  type(Pawtab_type),target,intent(in) :: Pawtab(Cryst%ntypat)
301  type(pawcprj_type),intent(in) :: in_Cprj(Cryst%natom,nspinor*nband_k)
302  type(pawcprj_type),intent(inout) :: out_Cprj(Cryst%natom,nspinor*nband_k) !vz_i
303 
304 !Local variables-------------------------------
305 !scalars
306  integer :: iatom,iat_sym,iband,ibsp_bz
307  integer :: ibsp_ibz,ik_ibz,indexj,ispinor,isym,itim
308  integer :: itypat,jl,jl0,jlmn,jln,jln0,jlpm,jm,jn,lmax,mm,ncpgr
309  real(dp) :: arg,wtk
310  logical :: isirred
311 !arrays
312  integer :: r0(3) !,nlmn_atom(Cryst%natom)
313  integer,ABI_CONTIGUOUS pointer :: indlmn(:,:)
314  real(dp) :: dum(2,nspinor),kbz(3),kirr(3),phase(2),swp(2),tmp(2,nspinor)
315  real(dp),allocatable :: DS_mmpl(:,:,:)
316 
317 ! *********************************************************************
318 
319  ncpgr = in_Cprj(1,1)%ncpgr
320  ABI_CHECK(ncpgr==0,"Derivatives of cprj are not coded")
321 
322 !Get the index of the IBZ image associated to the BZ k-point ik_bz and related simmetry.
323  call get_BZ_item(Kmesh,ik_bz,kbz,ik_ibz,isym,itim,isirred=isirred)
324 
325  if (isirred) then  ! It is a point in the IBZ, Symmetrization is not needed.
326    call pawcprj_copy(in_Cprj,out_Cprj)
327    RETURN  
328  end if
329 !
330 !The corresponding point kirr in the IBZ.
331  call get_IBZ_item(Kmesh,ik_ibz,kirr,wtk)
332 !
333 !=== DS_mmpl is the rotation matrix for real spherical harmonics associated to symrec(:,:,isym) ===
334 !* Note the convention used by Blanco in Eq. 27 : DS_mmp multiply spherical harmonics as row vectors
335  lmax=Pawang%l_max-1 ! l_max is Max l+1
336  ABI_ALLOCATE(DS_mmpl,(2*lmax+1,2*lmax+1,lmax+1))
337  DS_mmpl=Pawang%zarot(:,:,:,isym)
338 
339 !Local copy.
340 !do iatom=1,Cryst%natom
341 !nlmn_atom(iatom)=Pawtab(Cryst%typat(iatom))%lmn_size
342 !end do
343 !call pawcprj_alloc(out_Cprj,ncpgr,nlmn_atom)
344 !
345 !===========================================
346 !==== Loop over atoms to be symmetrized ====
347 !===========================================
348  do iatom=1,Cryst%natom
349    itypat=Cryst%typat(iatom)
350    iat_sym=Cryst%indsym(4,isym,iatom)
351    indlmn => Pawtab(itypat)%indlmn
352    r0=Cryst%indsym(1:3,isym,iatom) ! R^{-1} (xred(:,iatom)-tnons) = xred(:,iat_sym) + r0.
353    arg=two_pi*dot_product(kirr,r0)
354    phase(1)=COS(arg)
355    phase(2)=SIN(arg)
356 !
357 !  Loop over the (jl,jm,jn) components to be symmetrized.
358    jl0=-1; jln0=-1; indexj=1
359    do jlmn=1,Pawtab(itypat)%lmn_size
360      jl  =indlmn(1,jlmn)
361      jm  =indlmn(2,jlmn)
362      jn  =indlmn(3,jlmn)
363      jln =indlmn(5,jlmn)
364      jlpm=1+jl+jm
365      if (jln/=jln0) indexj=indexj+2*jl0+1
366 !
367 !    === For each band, calculate contribution due to rotated real spherical harmonics ===
368 !    FIXME check this expression; according to Blanco I should have D(S^-1} but it seems D(S) is correct
369 !    Recheck spinorial case, presently is wrong
370      ibsp_ibz=0
371      ibsp_bz=0
372      do iband=1,nband_k
373 
374        tmp(:,:)=zero
375        do ispinor=1,nspinor
376          ibsp_ibz=ibsp_ibz+1
377          do mm=1,2*jl+1
378            tmp(1,ispinor)=tmp(1,ispinor)+DS_mmpl(mm,jlpm,jl+1)*in_Cprj(iat_sym,ibsp_ibz)%cp(1,indexj+mm)
379            tmp(2,ispinor)=tmp(2,ispinor)+DS_mmpl(mm,jlpm,jl+1)*in_Cprj(iat_sym,ibsp_ibz)%cp(2,indexj+mm)
380          end do
381        end do !ispinor
382 !
383 !      * Apply the phase to account if the symmetric atom belongs to a different unit cell.
384        do ispinor=1,nspinor
385          dum(1,ispinor)=tmp(1,ispinor)*phase(1)-tmp(2,ispinor)*phase(2)
386          dum(2,ispinor)=tmp(1,ispinor)*phase(2)+tmp(2,ispinor)*phase(1)
387        end do
388 !
389 !      * If required, apply time-reversal symmetry to retrieve the correct point in the BZ.
390        if (itim==2) then
391          if (nspinor==1) then
392            dum(2,1)=-dum(2,1)
393          else if (nspinor==2) then ! TODO rotate wavefunction in spinor space.
394            swp(:)=dum(:,1)
395            dum(1,1)= dum(1,2)
396            dum(2,1)=-dum(2,2)
397            dum(1,2)=-swp(1)
398            dum(2,2)= swp(2)
399          end if
400        end if
401 !
402 !      ==== Save values ====
403        do ispinor=1,nspinor
404          ibsp_bz=ibsp_bz+1
405          out_Cprj(iatom,ibsp_bz)%cp(1,jlmn)=dum(1,ispinor)
406          out_Cprj(iatom,ibsp_bz)%cp(2,jlmn)=dum(2,ispinor)
407        end do
408      end do !iband
409 
410      jl0=jl; jln0=jln
411    end do !jlmn
412  end do !iatom
413 
414  ABI_DEALLOCATE(DS_mmpl)
415 
416 end subroutine paw_symcprj_op