TABLE OF CONTENTS


ABINIT/m_read_plowannier [ Modules ]

[ Top ] [ Modules ]

NAME

  m_read_plowannier

FUNCTION

  Read Wannier coefficient in the file forlb.ovlp for ucrpa calculation
  this file was typically created in a DFT run with usedmft=1 and nbandkss -1

COPYRIGHT

 Copyright (C) 2006-2018 ABINIT group (BAmadon)
 This file is distributed under the terms of the
 GNU General Public License, see ~abinit/COPYING
 or http://www.gnu.org/copyleft/gpl.txt .

INPUTS

OUTPUT

PARENTS

CHILDREN

SOURCE

26 #if defined HAVE_CONFIG_H
27 #include "config.h"
28 #endif
29 
30 #include "abi_common.h"
31 
32 MODULE m_read_plowannier
33 
34  use defs_basis
35  implicit none
36 
37  private
38 
39  public :: read_plowannier

m_read_plowannier/read_plowannier [ Modules ]

[ Top ] [ m_read_plowannier ] [ Modules ]

NAME

 read_plowannier

FUNCTION

  Read Wannier coefficient in the file forlb.ovlp for ucrpa calculation
  this file was typically created in a DFT run with usedmft=1 and nbandkss -1

COPYRIGHT

 Copyright (C) 2006-2018 ABINIT group (BAmadon)
 This file is distributed under the terms of the
 GNU General Public License, see ~abinit/COPYING
 or http://www.gnu.org/copyleft/gpl.txt .

INPUTS

 Cryst<cryst_t>= data type gathering info on symmetries and unit cell
    %natom=number of atoms
    %nsym=number of symmetries
    %xred(3,natom)=reduced coordinated of atoms
    %typat(natom)=type of each atom
    %rprimd(3,3)=dimensional primitive translations in real space (bohr)
    %timrev= 2 if time reversal can be used, 1 otherwise
 Kmesh <kmesh_t>= datatype gathering parameters related to the k-point sampling
    %nibz=number of k-points in the IBZ
    %nbz=number of k-points in the BZ
    %bz(3,nbz)=reduced coordinates for k-points in the full Brillouin zone
    %ibz(3,nibz)=reduced coordinates for k-points in the irreducible wedge
    %tab(nbz)=mapping between a kpt in the BZ (array bz) and the irred point in the array ibz
    %tabi(nbz)= -1 if inversion is needed to obtain this particular kpt in the BZ, 1 means identity
    %tabo(nbz)= for each point in the BZ, the index of the symmetry operation S in reciprocal
      space which rotates k_IBZ onto \pm k_BZ (depending on tabi)
    %tabp(nbz)= For each k_BZ, it gives the phase factors associated to non-symmorphic operations, i.e
      e^{-i 2 \pi k_IBZ \cdot R{^-1}t} == e{-i 2\pi k_BZ cdot t} where :
      \transpose R{-1}=S and (S k_IBZ) = \pm k_BZ (depending on ktabi)
    %tabr(nfftot,nbz) For each point r on the real mesh and for each k-point in the BZ, tabr
      gives the index of (R^-1 (r-t)) in the FFT array where R=\transpose S^{-1} and k_BZ=S k_IBZ.
      t is the fractional translation associated to R
  luwindow: T if ucrpa_window is activated, F if not.
  prtvol: integer to give the amount of printing.
 nsppol : number of spin polarization.
 Pawang<pawang_type> angular mesh discretization and related data:

OUTPUT

 bandinf, bandsup : lower and upper bands for define Wannier functions
 coeffW_BZ(nsppol,bandinf:bandsup,nvz,2*lcor+1)
 lcor : angular momentum for correlated orbitals
 itypatcor : correlated species

PARENTS

      cchi0,cchi0q0,prep_calc_ucrpa

CHILDREN

      get_bz_item,wrtout

SOURCE

101 subroutine read_plowannier(cryst,bandinf,bandsup,coeffW_BZ,itypatcor,Kmesh,lcor,luwindow,nspinor,nsppol,pawang,prtvol,ucrpa_bands)
102 
103  use defs_basis
104  use defs_datatypes
105  use defs_abitypes
106  use m_abicore
107  use m_errors
108 
109  use m_io_tools,      only : open_file
110  use m_crystal,       only : crystal_t
111  use m_bz_mesh,       only : kmesh_t, get_BZ_item
112  use m_pawang,        only : pawang_type
113 
114 !This section has been created automatically by the script Abilint (TD).
115 !Do not modify the following lines by hand.
116 #undef ABI_FUNC
117 #define ABI_FUNC 'read_plowannier'
118 !End of the abilint section
119 
120  implicit none
121 !Arguments ------------------------------------
122 !types and arrays
123  type(kmesh_t),intent(in) :: Kmesh
124  type(crystal_t),intent(in) :: Cryst
125  complex(dpc), allocatable, intent(inout) :: coeffW_BZ(:,:,:,:,:,:)
126  type(Pawang_type),intent(in) :: Pawang
127 !scalars
128  logical, intent(inout) :: luwindow
129  integer, intent(out) :: bandinf,bandsup,itypatcor,lcor
130  integer, intent(in) :: nspinor,nsppol,prtvol
131  integer, intent(in) :: ucrpa_bands(2)
132 
133 !Local variables-------------------------------
134  character(len=500) :: message,msg
135  integer :: at_indx,ik_ibz,band1,m1,m2,spin,ik_bz,dummy,isym,itim,iat,indx,ispinor,unt
136  real(dp) :: xx,yy
137  real(dp) :: kbz(3)
138  complex(dpc), allocatable :: coeffW_IBZ(:,:,:,:,:,:)
139 ! *********************************************************************
140  write(message,*) "Read wannier in iBZ"
141  call wrtout(std_out,message,'COLL')
142 
143  if (open_file('forlb.ovlp',msg,newunit=unt,form='formatted',status='unknown') /= 0) then
144    MSG_ERROR(msg)
145  end if
146  rewind(unt)
147  read(unt,*) message, lcor,itypatcor
148  read(unt,*) message, bandinf,bandsup
149  write(std_out,*) 'read from forlb.ovlp',lcor, bandinf,bandsup
150  write(std_out,*) "for wannier", bandinf,bandsup
151  if(prtvol>0) then
152  endif
153 
154  if(.not.luwindow.and.(ucrpa_bands(1)/=bandinf.or.ucrpa_bands(2)/=bandsup)) then
155    write(msg,'(a,a)')' Bands used for Wannier construction and cRPA construction differ',&
156 & 'It might be physically correct and it is possible with the current implementation'
157    call wrtout(std_out,msg,'COLL')
158    call wrtout(ab_out,msg,'COLL')
159 !   MSG_ERROR(msg)
160  end if
161 
162 !Do not dead the bandinf, bandinf redondance information
163 
164  ABI_ALLOCATE(coeffW_IBZ,(Cryst%nattyp(itypatcor),nsppol,bandinf:bandsup,Kmesh%nibz,nspinor,2*lcor+1))
165  coeffW_IBZ=czero
166  do spin=1,nsppol
167    do ik_ibz=1,Kmesh%nibz
168      !read k
169       read(unt,*)
170       do band1=bandinf,bandsup
171      !read band
172         read(unt,*)
173      !read projection
174         do ispinor=1,nspinor
175           do iat=1,Cryst%nattyp(itypatcor)
176             do m1=1,2*lcor+1
177               read(unt,*) dummy,dummy,dummy,xx,yy
178               coeffW_IBZ(iat,spin,band1,ik_ibz,ispinor,m1)=cmplx(xx,yy)
179             end do
180           end do
181         end do
182       end do
183    end do
184  end do
185  close(unt)
186 
187  ABI_ALLOCATE(coeffW_BZ,(Cryst%nattyp(itypatcor),nsppol,bandinf:bandsup,Kmesh%nbz,nspinor,2*lcor+1))
188  coeffW_BZ=czero
189 
190  if (Kmesh%nbz==Kmesh%nibz) then
191    coeffW_BZ=coeffW_IBZ
192  else if (Cryst%nsym==1) then
193    write(message,*) "Reconstruct in full BZ"
194    call wrtout(std_out,message,'COLL')
195    call wrtout(ab_out,message,'COLL')
196    do ik_bz=1,Kmesh%nbz
197 !     if(prtvol>=10) write(6,*) "ik",ik_bz,Kmesh%tab(ik_bz),Kmesh%tabi(ik_bz),Kmesh%tabo(ik_bz)
198      if(Kmesh%tabi(ik_bz)==1) then
199        coeffW_BZ(:,:,:,ik_bz,:,:)=coeffW_IBZ(:,:,:,Kmesh%tab(ik_bz),:,:)
200      else if(Kmesh%tabi(ik_bz)==-1) then
201        coeffW_BZ(:,:,:,ik_bz,:,:)=conjg(coeffW_IBZ(:,:,:,Kmesh%tab(ik_bz),:,:))
202      endif
203 !     if(prtvol>=10)write(6,*) "coeffW 21",coeffW_BZ(1,1,bandinf,ik_bz,:)
204 !     if(prtvol>=10)write(6,*) "coeffW 21",coeffW_BZ(1,1,bandinf+1,ik_bz,:)
205 !     if(prtvol>=10)write(6,*) "coeffW 25",coeffW_BZ(1,1,bandsup,ik_bz,:)
206    enddo
207  else if (Cryst%nsym>1) then
208    write(message,*) "Reconstruct in full BZ (nsym=0)"
209    call wrtout(std_out,message,'COLL')
210    do ik_bz=1,Kmesh%nbz
211    !write(6,*) ik_bz,Kmesh%nbz
212      call get_BZ_item(Kmesh,ik_bz,kbz,ik_ibz,isym,itim)
213      do indx=1,cryst%nattyp(itypatcor)
214        iat=cryst%atindx1(indx) ! correct index for the full atom list
215        at_indx=cryst%indsym(4,isym,iat) !! see eg sym_matlu and m_crystal
216        do spin=1,nsppol
217          do ispinor=1,nspinor
218            do m1=1,2*lcor+1
219              do m2=1,2*lcor+1
220              coeffW_BZ(indx,spin,:,ik_bz,ispinor,m1)= coeffW_BZ(indx,spin,:,ik_bz,ispinor,m1)&
221 &                     +coeffW_IBZ(at_indx,spin,:,ik_ibz,ispinor,m2)*pawang%zarot(m2,m1,lcor+1,isym)
222              enddo
223 !             write(6,'(20f7.3)') (pawang%zarot(m1,m2,lcor+1,isym),m2=1,2*lcor+1)
224            enddo
225          enddo
226        enddo
227      enddo
228 !     do m1=1,2*lcor+1
229 !      write(6,*) "coeffW_IBZ",coeffW_IBZ(1,8,ik_ibz,m1)
230 !     enddo
231 !     do m1=1,2*lcor+1
232 !      write(6,*) "coeffW_ BZ",coeffW_BZ(1,8,ik_bz,m1)
233 !     enddo
234    enddo
235  end if
236  ABI_DEALLOCATE(coeffW_IBZ)
237 
238 
239 end subroutine read_plowannier