TABLE OF CONTENTS


ABINIT/mkfskgrid [ Functions ]

[ Top ] [ Functions ]

NAME

 mkfskgrid

FUNCTION

 This routine sets up the full FS kpt grid by symmetry

COPYRIGHT

 Copyright (C) 2004-2018 ABINIT group (MVer)
 This file is distributed under the terms of the
 GNU General Public Licence, see ~abinit/COPYING
 or http://www.gnu.org/copyleft/gpl.txt .
 For the initials of contributors, see ~abinit/doc/developers/contributors.txt .

INPUTS

  nsym    = number of symmetries for the full system
  symrec  = reciprocal space symmetries (those for the kpts)
  timrev  = 1 if time reversal symmetry is to be used

OUTPUT

  elph_k datastructure:
  elph_k%nkpt           = full number of kpoints close to the FS
  elph_k%kpt            = full set of kpoints close to the FS
  elph_k%wtkirr         = weights of the irreducible kpoints
  elph_k%kphon_irr2full = indices of irred kpoints in full array

NOTES

  WARNING: supposes kpt grid has full symmetry!! Not always true!!!
    but should be for Monkhorst-Pack, efficient grids.
    otherwise you get an error message in interpolate_gkk because
    an FS kpt can not be found in the gkk file.

PARENTS

      elphon

CHILDREN

      destroy_kptrank,get_rank_1kpt,mkkptrank,sort_int,wrap2_pmhalf,wrtout

SOURCE

 43 #if defined HAVE_CONFIG_H
 44 #include "config.h"
 45 #endif
 46 
 47 #include "abi_common.h"
 48 
 49 subroutine mkFSkgrid (elph_k, nsym, symrec, timrev)
 50 
 51  use defs_basis
 52  use defs_elphon
 53  use m_kptrank
 54  use m_profiling_abi
 55  use m_errors
 56  use m_sort
 57 
 58  use m_numeric_tools,   only : wrap2_pmhalf
 59 
 60 !This section has been created automatically by the script Abilint (TD).
 61 !Do not modify the following lines by hand.
 62 #undef ABI_FUNC
 63 #define ABI_FUNC 'mkFSkgrid'
 64  use interfaces_14_hidewrite
 65 !End of the abilint section
 66 
 67  implicit none
 68 
 69 !Arguments ------------------------------------
 70 !scalars
 71  integer,intent(in) :: nsym,timrev
 72  type(elph_kgrid_type),intent(inout) :: elph_k
 73 !arrays
 74  integer,intent(in) :: symrec(3,3,nsym)
 75 
 76 !Local variables-------------------------------
 77 !scalars
 78  integer :: ikpt1,ikpt2,isym,itim,new,symrankkpt
 79  real(dp) :: timsign, res
 80  character(len=500) :: message
 81 
 82 !arrays
 83  real(dp) :: kpt(3),redkpt(3)
 84  integer, allocatable :: sortindexing(:), rankallk(:)
 85 
 86  integer, allocatable :: tmpkphon_full2irr(:,:)
 87  real(dp), allocatable :: tmpkpt(:,:)
 88 
 89 ! *************************************************************************
 90 
 91  if(timrev /= 1 .and. timrev /= 0)then
 92    write (message,'(a,i0)')' timrev must be 1 or 0 but found timrev= ',timrev
 93    MSG_BUG(message)
 94  end if
 95 
 96  ABI_ALLOCATE(tmpkphon_full2irr,(3,2*elph_k%nkptirr*nsym))
 97  tmpkphon_full2irr = -1
 98 
 99  ABI_ALLOCATE(tmpkpt,(3,2*elph_k%nkptirr*nsym))
100 
101  ABI_ALLOCATE(elph_k%wtkirr,(elph_k%nkptirr))
102  elph_k%wtkirr(:) = zero
103 
104 !first allocation for irred kpoints - will be destroyed below
105  call mkkptrank (elph_k%kptirr,elph_k%nkptirr,elph_k%kptrank_t)
106  ABI_ALLOCATE(rankallk,(elph_k%kptrank_t%max_rank))
107 
108 !elph_k%kptrank_t%invrank is used as a placeholder in the following loop
109  rankallk = -1
110  elph_k%kptrank_t%invrank = -1
111 
112 !replicate all irred kpts by symmetry to get the full k grid.
113  elph_k%nkpt=0 !zero k-points found so far
114  do isym=1,nsym
115    do itim=0,1
116      timsign = one-two*itim
117      do ikpt1=1,elph_k%nkptirr
118 !      generate symmetrics of kpt ikpt1
119        kpt(:) = timsign*(symrec(:,1,isym)*elph_k%kptirr(1,ikpt1) + &
120 &       symrec(:,2,isym)*elph_k%kptirr(2,ikpt1) + &
121 &       symrec(:,3,isym)*elph_k%kptirr(3,ikpt1))
122        
123        call get_rank_1kpt (kpt,symrankkpt,elph_k%kptrank_t)
124 
125 !      is the kpt on the full grid (may have lower symmetry than full spgroup)
126 !      is kpt among the full FS kpts found already?
127        if (elph_k%kptrank_t%invrank(symrankkpt) == -1) then
128          elph_k%wtkirr(ikpt1)=elph_k%wtkirr(ikpt1)+1
129          elph_k%nkpt=elph_k%nkpt+1
130 
131          call wrap2_pmhalf(kpt(1),redkpt(1),res)
132          call wrap2_pmhalf(kpt(2),redkpt(2),res)
133          call wrap2_pmhalf(kpt(3),redkpt(3),res)
134          tmpkpt(:,elph_k%nkpt) = redkpt
135          tmpkphon_full2irr(1,elph_k%nkpt) = ikpt1
136 !        save sym that sends irred kpt ikpt1 onto full kpt
137          tmpkphon_full2irr(2,elph_k%nkpt) = isym
138          tmpkphon_full2irr(3,elph_k%nkpt) = itim
139          
140          elph_k%kptrank_t%invrank(symrankkpt) = elph_k%nkpt
141          rankallk(elph_k%nkpt) = symrankkpt
142        end if
143 
144      end do !end loop over irred k points
145    end do !end loop over timrev
146  end do !end loop over symmetry
147 
148  write(message,'(a,i0)')'mkfskgrid: after first evaluation, elph_k%nkpt= ', elph_k%nkpt
149  call wrtout(std_out,message,"COLL")
150 
151  elph_k%wtkirr(:) = elph_k%wtkirr(:) / elph_k%nkpt
152 
153 !copy the kpoints and full --> irred kpt map
154 !reorder the kpts to get rank increasing monotonically with a sort
155 !also reorder tmpkphon_full2irr
156  ABI_ALLOCATE(elph_k%kpt,(3,elph_k%nkpt))
157  ABI_ALLOCATE(elph_k%full2irr,(3,elph_k%nkpt))
158  ABI_ALLOCATE(sortindexing,(elph_k%nkpt))
159 
160  do ikpt1=1,elph_k%nkpt
161    sortindexing(ikpt1)=ikpt1
162  end do
163  call sort_int(elph_k%nkpt, rankallk, sortindexing)
164  do ikpt1=1,elph_k%nkpt
165    if (sortindexing(ikpt1) < 1 .or. sortindexing(ikpt1) > elph_k%nkpt) then
166      MSG_BUG('sorted k ranks are out of bounds: 1 to nkpt')
167    end if 
168    elph_k%kpt(:,ikpt1) = tmpkpt(:,sortindexing(ikpt1))
169    elph_k%full2irr(:,ikpt1) = tmpkphon_full2irr(:,sortindexing(ikpt1))
170  end do
171 
172  ABI_DEALLOCATE(sortindexing)
173  ABI_DEALLOCATE(rankallk)
174  ABI_DEALLOCATE(tmpkphon_full2irr)
175  ABI_DEALLOCATE(tmpkpt)
176  call destroy_kptrank (elph_k%kptrank_t)
177 
178 
179 !make proper full rank arrays
180  call mkkptrank (elph_k%kpt,elph_k%nkpt,elph_k%kptrank_t)
181 
182 
183 !find correspondence table between irred FS kpoints and a full one
184  ABI_ALLOCATE(elph_k%irr2full,(elph_k%nkptirr))
185  elph_k%irr2full(:) = 0
186 
187  do ikpt1=1,elph_k%nkptirr
188    call get_rank_1kpt (elph_k%kptirr(:,ikpt1),symrankkpt,elph_k%kptrank_t)
189    elph_k%irr2full(ikpt1) = elph_k%kptrank_t%invrank(symrankkpt)
190  end do
191 
192 !find correspondence table between FS kpoints under symmetry
193  ABI_ALLOCATE(elph_k%full2full,(2,nsym,elph_k%nkpt))
194  elph_k%full2full(:,:,:) = -999
195 
196  do ikpt1=1,elph_k%nkpt
197 !  generate symmetrics of kpt ikpt1
198    do isym=1,nsym
199      do itim=0,timrev
200        timsign = one-two*itim
201        kpt(:) = timsign*(symrec(:,1,isym)*elph_k%kpt(1,ikpt1) + &
202 &       symrec(:,2,isym)*elph_k%kpt(2,ikpt1) + &
203 &       symrec(:,3,isym)*elph_k%kpt(3,ikpt1))
204 
205 !      which kpt is it among the full FS kpts
206        call get_rank_1kpt (kpt,symrankkpt,elph_k%kptrank_t)
207        ikpt2 = elph_k%kptrank_t%invrank(symrankkpt)
208        new=1
209        if (ikpt2 /= -1) then
210          elph_k%full2full(itim+1,isym,ikpt2) = ikpt1
211          new = 0
212        end if
213 
214        if (new == 1) then
215          write(std_out,*) ' mkfskgrid Error: FS kpt ',ikpt1,' has no symmetric under sym', isym,' with itim ',itim
216          write(std_out,*) ' redkpt = ', redkpt
217          write(std_out,*) ' symrankkpt,ikpt2 = ', symrankkpt,ikpt2
218          MSG_ERROR("Fatal error, cannot continue")
219        end if
220      end do
221    end do
222  end do
223 
224 !got nkpt, tmpkpt, kphon_full2irr, kphon_full2full, and wtkirr
225 
226 end subroutine mkFSkgrid