TABLE OF CONTENTS


ABINIT/mlwfovlp_setup [ Functions ]

[ Top ] [ Functions ]

NAME

 mlwfovlp_setup

FUNCTION

 Routine which creates table g1 and ovikp  necessary to compute
 overlap for Wannier code (www.wannier.org f90 version).

COPYRIGHT

 Copyright (C) 2005-2018 ABINIT group (BAmadon,FJollet)
 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

  atom_symbols(natom)= table of symbol for each atom
                                          and each |p_lmn> non-local projector
  dtset <type(dataset_type)>=all input variables for this dataset
  filew90_win(nsppol) secondary input files for w90
  lwanniersetup= flag: only 1 is fully working.
  natom              =number of atoms in cell.
  mband=maximum number of bands
  natom=number of atoms in cell.
  nkpt=number of k points.
  num_bands(isppol)=number of bands actually used to construct the wannier function
  nwan(isppol)= number of wannier fonctions (read in wannier90.win).
  dtset <type(dataset_type)>=all input variables for this dataset
  real_lattice(3,3)=dimensional primitive translations for real space
                 in format required by wannier90
  recip_lattice(3,3)=dimensional primitive translations for reciprocal space
                 in format required by wannier90
  rprimd(3,3)=dimensional primitive translations for real space (bohr)
  seed_name=character string for generating wannier90 filenames
  spin = just used for nsppol>1 ; 0 both, 1 just spin up, 2 just spin down
  xcart(3,natom)=atomic coordinates in bohr
  xred(3,natom)=reduced dimensionless atomic coordinates

OUTPUT

  band_in(mband,nsppol)   = band to take into account for wannier calculation
  g1(3,nkpt,nntot) = G vector shift which is necessary to obtain k1+b
                     from k2 in the case where k1+b does not belong to the 1st BZ.
  nband_inc(nsppol) = # of included bands
  nntot            = number of k-point neighbour
  ovikp(nkpt,nntot)= gives  nntot value of k2 (in the BZ) for each k1  (k2=k1+b mod(G))

SIDE EFFECTS

  (only writing, printing)

NOTES

PARENTS

      mlwfovlp

CHILDREN

      atomdata_from_znucl,wannier_setup,wrtout

SOURCE

 60 #if defined HAVE_CONFIG_H
 61 #include "config.h"
 62 #endif
 63 
 64 #include "abi_common.h"
 65 
 66 
 67  subroutine mlwfovlp_setup(atom_symbols,band_in,dtset,filew90_win,gamma_only,&
 68 & g1,lwanniersetup,mband,natom,nband_inc,nkpt,&
 69 & nntot,num_bands,num_nnmax,nsppol,nwan,ovikp,&
 70 & proj_l,proj_m,proj_radial,proj_site,proj_x,proj_z,proj_zona,&
 71 & real_lattice,recip_lattice,rprimd,seed_name,spin,spinors,xcart,xred)
 72 
 73  use defs_basis
 74  use defs_abitypes
 75  use defs_wannier90
 76  use m_errors
 77  use m_profiling_abi
 78  use m_atomdata
 79 
 80  use m_io_tools,   only : open_file
 81 
 82 !This section has been created automatically by the script Abilint (TD).
 83 !Do not modify the following lines by hand.
 84 #undef ABI_FUNC
 85 #define ABI_FUNC 'mlwfovlp_setup'
 86  use interfaces_14_hidewrite
 87 !End of the abilint section
 88 
 89  implicit none
 90 
 91 !Arguments---------------------------
 92 ! scalars
 93 !scalars
 94  integer,intent(in) :: lwanniersetup,mband,natom,nkpt,nsppol
 95  integer,intent(in) :: num_nnmax,spin
 96  integer,intent(out) :: nband_inc(nsppol),nntot,num_bands(nsppol),nwan(nsppol)
 97  logical,intent(in) :: gamma_only,spinors
 98  type(dataset_type),intent(in) :: dtset
 99 !arrays
100  integer,intent(out) :: g1(3,nkpt,num_nnmax),ovikp(nkpt,num_nnmax)
101  integer,intent(out) :: proj_l(mband,nsppol),proj_m(mband,nsppol),proj_radial(mband,nsppol)
102  real(dp),intent(in) :: real_lattice(3,3)
103  real(dp),intent(in) :: recip_lattice(3,3),rprimd(3,3),xred(3,natom)
104  real(dp),intent(out) :: proj_site(3,mband,nsppol),proj_x(3,mband,nsppol),proj_z(3,mband,nsppol)
105  real(dp),intent(out) :: proj_zona(mband,nsppol),xcart(3,natom)
106  logical,intent(out) :: band_in(mband,nsppol)
107  character(len=3),intent(out) :: atom_symbols(natom)
108  character(len=fnlen),intent(in) :: seed_name(nsppol),filew90_win(nsppol)
109 
110 !Local variables---------------------------
111 !scalars
112  integer :: iatom,icb,ikpt,ikpt1,intot,isppol,itypat,jj,mband_,unt
113  real(dp) :: znucl1
114  character(len=2) :: symbol
115  character(len=500) :: message
116  character(len=fnlen) :: filew90_nnkp
117  type(atomdata_t) :: atom
118 !arrays
119  integer :: exclude_bands(mband,nsppol),ngkpt(3)
120 
121 ! *************************************************************************
122 
123 !^^^^^^^^^^^^^^^^read wannier90.nnkp^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
124  if(lwanniersetup==0) then  !this part is not coded for nsppol>1
125    isppol=1
126    filew90_nnkp=trim(seed_name(isppol))//'.nnkp'
127    if (open_file(filew90_nnkp,message,newunit=unt,form='formatted',status='old') /= 0) then
128      MSG_ERROR(message)
129    end if
130    read(unt,*)
131    read(unt,*) nntot , mband_, nwan(1)
132    write(message, '(a,a,i6,i6,i6)' )ch10,&
133 &   ' mlwfovlp_setup nntot,mband,nwan ', nntot,mband_,nwan(1)
134    call wrtout(std_out,message,'COLL')
135    if(mband_.ne.mband) then
136      write(message, '(4a)' )&
137 &     'mband_ is not equal to mband ',ch10,&
138 &     'Action: check ',trim(filew90_nnkp)
139      MSG_ERROR(message)
140    end if
141    if(nwan(1)>mband) then
142      write(message, '(4a)' )&
143 &     'nwan > mband ',ch10,&
144 &     'Action: check ',trim(filew90_nnkp)
145      MSG_ERROR(message)
146    end if
147    if(nwan(1)==0) then
148      write(message, '(4a)' )&
149 &     'nwan = 0 ',ch10,&
150 &     'Action: check ',trim(filew90_nnkp)
151      MSG_ERROR(message)
152    end if
153    do ikpt=1,nkpt
154      do intot=1,nntot
155 !      ikpt1: k point  (ikpt=ikpt1)
156 !      ovikp(intot,ikpt): neighbour number intot for ikpt
157 !      g1(1:3,intot,ikpt): non reciprocal space vector between the 2 k-points
158        read(unt,*)  &
159 &       ikpt1,ovikp(ikpt,intot),(g1(jj,ikpt,intot),jj=1,3)
160        if(ikpt1.ne.ikpt) write(std_out,*) "warning: ikpt1 .ne ikpt : ?"
161      end do
162    end do
163    close(unt)
164    write(message, '(3a)' )ch10,&
165 &   trim(filew90_nnkp),'wannier90.nnkp has been read !'
166    call wrtout(std_out,message,'COLL')
167 
168    message = ' exclude bands is not given in this case (not implemented) '
169    MSG_ERROR(message)
170 
171 !  ^^^^^^^^^^^^^^^^^^^^^^^ call wannier_setup begin^^^^^^^^^^^^^^^^^^^^^^^^
172  else if (lwanniersetup==1) then
173    num_bands(:)=mband
174 !  num_nnmax=12 !limit fixed for compact structure in wannier_setup.
175    ovikp=0.d0
176 !  "When nshiftk=1, kptrlatt is initialized as a diagonal (3x3) matrix, whose diagonal 
177 !  elements are the three values ngkpt(1:3)"
178    ngkpt(1)=dtset%kptrlatt(1,1)
179    ngkpt(2)=dtset%kptrlatt(2,2) !  have to verif kptrlatt is diagonal
180    ngkpt(3)=dtset%kptrlatt(3,3)
181    do iatom=1,natom
182      itypat=dtset%typat(iatom)
183      znucl1=dtset%znucl(itypat)
184      call atomdata_from_znucl(atom, znucl1) 
185      symbol=trim(adjustl(atom%symbol))
186 !    write(309,*) symbol
187      atom_symbols(iatom)=symbol
188      xcart(:,iatom)=rprimd(:,1)*xred(1,iatom)+&
189 &     rprimd(:,2)*xred(2,iatom)+&
190 &     rprimd(:,3)*xred(3,iatom)
191    end do ! iatom
192 !  write(std_out,*) xcart
193 !  write(std_out,*) Bohr_Ang
194 !  write(std_out,*) rprimd*Bohr_Ang
195 !  write(std_out,*) seed_name
196 !  write(std_out,*) ngkpt
197 !  write(std_out,*) nkpt
198 !  write(std_out,*) mband
199 !  write(std_out,*) natom
200 !  write(std_out,*) atom_symbols
201    write(message, '(a,a)' )ch10,&
202 &   '** mlwfovlp_setup:  call wannier90 library subroutine wannier_setup'
203    call wrtout(std_out,message,'COLL')
204 #if defined HAVE_WANNIER90
205    nwan(:)=0
206    num_bands(:)=0
207    do isppol=1,nsppol
208      if(spin.ne.0 .and. spin.ne.isppol) cycle
209      call wannier_setup(seed_name(isppol),ngkpt,nkpt&                    !input
210 &    ,real_lattice,recip_lattice,dtset%kpt&                      !input
211 &    ,mband,natom,atom_symbols,xcart*Bohr_Ang&                   !input
212 &    ,gamma_only,spinors&                                        !input
213 &    ,nntot,ovikp,g1,num_bands(isppol),nwan(isppol)&                     !output
214 &    ,proj_site(:,:,isppol),proj_l(:,isppol)&                    !output
215 &    ,proj_m(:,isppol),proj_radial(:,isppol)&                    !output
216 &    ,proj_z(:,:,isppol),proj_x(:,:,isppol)&                     !output
217 &    ,proj_zona(:,isppol),exclude_bands(:,isppol) )               !output
218    end do !isppol
219 #else
220    ABI_UNUSED(gamma_only)
221    ABI_UNUSED(real_lattice)
222    ABI_UNUSED(recip_lattice)
223    ABI_UNUSED(spinors)
224 #endif
225 !  do isppol=1,nsppol
226 !  if(spin.ne.0 .and. spin.ne.isppol) cycle
227 !  write(std_out,*)  "1", nntot,nwan(isppol)
228 !  write(std_out,*)  "2",num_bands(isppol)  ! states on which wannier functions are computed
229 !  write(std_out,*)  "3", proj_site(:,1:nwan(isppol),isppol)
230 !  write(std_out,*)  "4",proj_l(1:nwan(isppol),isppol)
231 !  write(std_out,*)  "5",proj_m(1:nwan(isppol),isppol)
232 !  write(std_out,*)  "6", proj_radial(1:nwan(isppol),isppol)
233 !  write(std_out,*)  "7", proj_z(:,1:nwan(isppol),isppol)
234 !  write(std_out,*)  "8", proj_x(:,1:nwan(isppol),isppol)
235 !  write(std_out,*)  "9",proj_zona(1:nwan(isppol),isppol)
236 !  write(std_out,*)  "10",exclude_bands(:,isppol)
237 !  end do!isppol
238 !  testdebug:  ovikp(1,1)=1
239 !  ^^^^^^^^^^^^^^^^^^^^^^^ end ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
240  end if  ! lwanniersetup
241  do isppol=1,nsppol
242    if(spin.ne.0 .and. spin.ne.isppol) cycle
243    band_in(:,isppol)=.true.
244    do icb=1,mband
245      if(exclude_bands(icb,isppol).ne.0)  band_in(exclude_bands(icb,isppol),isppol)=.false. 
246    end do
247    nband_inc(isppol)=0
248    do icb=1, mband
249      if (band_in(icb,isppol)) then
250        nband_inc(isppol)=nband_inc(isppol)+1
251      end if
252    end do
253  end do !isppol
254  if(any(mband.gt.num_bands(:))) then
255    write(message, '(a,a)' )ch10,&
256 &   '   Following bands are excluded from the calculation of wannier functions:'
257    call wrtout(std_out,message,'COLL')
258    
259    do isppol=1,nsppol
260      if(spin.ne.0 .and. spin.ne.isppol) cycle
261      if(nsppol==2) then
262        write(message,'("For spin",i4)')isppol
263 !      write(message,'(a,i)')'For spin=',isppol
264        call wrtout(std_out,message,'COLL')
265      end if !nsppol
266      do jj=1,mband-num_bands(isppol),10
267        write(message,'(10i7)') exclude_bands(jj:min(jj+9,mband-num_bands(isppol)),isppol)
268        call wrtout(std_out,message,'COLL')
269      end do
270    end do !isppol
271  end if
272 
273  do isppol=1,nsppol
274    if(spin.ne.0 .and. spin.ne.isppol) cycle
275    if(nsppol==2) then
276      write(message,'("For spin",i4)')isppol
277      call wrtout(std_out,message,'COLL')
278    end if !nsppol
279    write(message, '(a,i6,3a)' )ch10,&
280 &   nwan(isppol),' wannier functions will be computed (see ',trim(filew90_win(isppol)),')'
281    call wrtout(std_out,message,'COLL') 
282 !  write(std_out,*) exclude_bands(icb),band_in(icb)
283 !  ^^^^^^^^^^^^^^^END OF READING
284    write(message, '(a,i6,a)' )ch10,&
285 &   num_bands(isppol),' bands will be used to extract wannier functions'
286    call wrtout(std_out,message,'COLL')
287    if(num_bands(isppol).lt.nwan(isppol)) then
288      write(message, '(4a)' )&
289 &     ' number of bands is lower than the number of wannier functions',ch10,&
290 &     ' Action : check input file and ',trim(filew90_win(isppol))
291      MSG_ERROR(message)
292    else if (num_bands(isppol)==nwan(isppol)) then
293      write(message, '(a,a,a,a)' )ch10,&
294 &     '   Number of bands is equal to the number of wannier functions',ch10,&
295 &     '   Disentanglement will not be necessary'
296      call wrtout(std_out,message,'COLL')
297    else if  (num_bands(isppol).gt.nwan(isppol)) then
298      write(message, '(a,a,a,a)' )ch10,&
299 &     '   Number of bands is larger than the number of wannier functions',ch10,&
300 &     '   Disentanglement will be necessary'
301      call wrtout(std_out,message,'COLL')
302    end if
303    write(message, '(2x,a,a,i3,1x,a)' )ch10,&
304 &   '   Each k-point has', nntot,'neighbours'
305    call wrtout(std_out,message,'COLL')
306 
307  end do !isppol
308 
309 end subroutine mlwfovlp_setup