TABLE OF CONTENTS


ABINIT/m_fstab [ Modules ]

[ Top ] [ Modules ]

NAME

FUNCTION

  Tools for the management of a set of Fermi surface k-points

COPYRIGHT

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

PARENTS

SOURCE

18 #if defined HAVE_CONFIG_H
19 #include "config.h"
20 #endif
21 
22 #include "abi_common.h"
23 
24 
25 module m_fstab
26 
27  use defs_basis
28  use defs_abitypes
29  use m_abicore
30  use m_xmpi
31  use m_errors
32  use m_kptrank
33  use m_tetrahedron
34  use m_ebands
35 
36  use m_time,           only : cwtime
37  use m_fstrings,       only : itoa, sjoin
38  use m_numeric_tools,  only : bisect
39  use m_symtk,          only : matr3inv
40  use defs_datatypes,   only : ebands_t
41  use m_crystal,        only : crystal_t
42  use m_special_funcs,  only : dirac_delta
43  use m_kpts,           only : kpts_timrev_from_kptopt, listkk, smpbz
44 
45  implicit none
46 
47  private

m_fstab/fstab_findkg0 [ Functions ]

[ Top ] [ m_fstab ] [ Functions ]

NAME

  fstab_findkg0

FUNCTION

  Return the index `ikfs` of the k-point `kpt` in the FS-BZ. Return -1 if not found.

INPUTS

  kpt(3)=K-point in reduced coordinates

OUTPUT

   g0=Reciprocal lattice vector such that kpt = fstab%kpts(:, ikfs) + g0

PARENTS

CHILDREN

SOURCE

545 integer function fstab_findkg0(fstab, kpt, g0) result(ikfs)
546 
547 
548 !This section has been created automatically by the script Abilint (TD).
549 !Do not modify the following lines by hand.
550 #undef ABI_FUNC
551 #define ABI_FUNC 'fstab_findkg0'
552 !End of the abilint section
553 
554  implicit none
555 
556 !Arguments ------------------------------------
557 !scalars
558  type(fstab_t),intent(in) :: fstab
559 !arrays
560  integer,intent(out) :: g0(3)
561  real(dp),intent(in) :: kpt(3)
562 
563 ! *************************************************************************
564 
565  ikfs = kptrank_index(fstab%krank, kpt)
566  if (ikfs /= -1) then
567    g0 = nint(kpt - fstab%kpts(:, ikfs))
568  else
569    g0 = huge(1)
570  end if
571 
572 end function fstab_findkg0

m_fstab/fstab_free [ Functions ]

[ Top ] [ m_fstab ] [ Functions ]

NAME

  fstab_free

FUNCTION

  Free memory

INPUTS

OUTPUT

PARENTS

      m_phgamma

CHILDREN

      wrtout

SOURCE

165 subroutine fstab_free(fstab)
166 
167 
168 !This section has been created automatically by the script Abilint (TD).
169 !Do not modify the following lines by hand.
170 #undef ABI_FUNC
171 #define ABI_FUNC 'fstab_free'
172 !End of the abilint section
173 
174  implicit none
175 
176 !Arguments ------------------------------------
177  type(fstab_t),intent(inout) :: fstab
178 
179 ! ************************************************************************
180 
181  !@fstab_t
182 
183  ! integer
184  if (allocated(fstab%istg0)) then
185    ABI_FREE(fstab%istg0)
186  end if
187  if (allocated(fstab%bstcnt_ibz)) then
188    ABI_FREE(fstab%bstcnt_ibz)
189  end if
190 
191  ! real
192  if (allocated(fstab%kpts)) then
193    ABI_FREE(fstab%kpts)
194  end if
195  if (allocated(fstab%tetra_wtk)) then
196    ABI_FREE(fstab%tetra_wtk)
197  end if
198  if (allocated(fstab%tetra_wtk_ene)) then
199    ABI_FREE(fstab%tetra_wtk_ene)
200  end if
201 
202  ! types
203  call destroy_kptrank(fstab%krank)
204 
205 end subroutine fstab_free

m_fstab/fstab_init [ Functions ]

[ Top ] [ m_fstab ] [ Functions ]

NAME

  fstab_init

FUNCTION

  Initialize the tables for the FS integration.

INPUTS

  ebands<ebands_t>=The object describing the band structure.
  cryst<crystal_t>=Info on the crystalline structure.
  fsewin=Energy window in Hartree. Only states in [efermi-fsewin, efermi+fsewin] are included.
  integ_method=Flag selecting the integration method.
  kptrlatt(3,3)=k-point lattice specification
  nshiftk= number of shift vectors.
  shiftk(3,nshiftk)=shift vectors for k point generation
  comm=MPI communicator.

OUTPUT

  fstab(nsppol)=Tables with the correspondence between points of the Fermi surface (FS)
     and the k-points in ebands_t.

TODO

  Use a different algorithm to select k-points if tetra. First compute tetra weights
  then k-points contributing to FS integral are selected according to some threshold.

PARENTS

      m_phgamma

CHILDREN

      wrtout

SOURCE

243 subroutine fstab_init(fstab, ebands, cryst, fsewin, integ_method, kptrlatt, nshiftk, shiftk, comm)
244 
245 
246 !This section has been created automatically by the script Abilint (TD).
247 !Do not modify the following lines by hand.
248 #undef ABI_FUNC
249 #define ABI_FUNC 'fstab_init'
250 !End of the abilint section
251 
252  implicit none
253 
254 !Arguments ------------------------------------
255 !scalars
256  integer,intent(in) :: nshiftk,integ_method,comm
257  real(dp),intent(in) :: fsewin
258  type(ebands_t),intent(in) :: ebands
259  type(crystal_t),intent(in) :: cryst
260 !arrays
261  integer,intent(in) :: kptrlatt(3,3)
262  real(dp),intent(in) :: shiftk(3,nshiftk)
263  type(fstab_t),target,intent(out) :: fstab(ebands%nsppol)
264 
265 !Local variables-------------------------------
266 !scalars
267  integer,parameter :: option0=0,brav1=1,bcorr0=0
268  integer :: nkfs,spin,band,nband_k,i1,i2,ib,blow,ik_full,ik_ibz,nkibz,sppoldbl,timrev
269  integer :: ik,mkpt,nkpt_full,ierr !ikfull,my_kstart,my_kstop,isym,itime,tsign
270  integer :: bstart_k,bstop_k,nene,ifermi,bmin,bmax
271  real(dp) :: elow,ehigh,ebis,enemin,enemax,deltaene,max_occ,dksqmax,cpu,wall,gflops
272  logical :: inwin
273  character (len=80) :: errstr
274 !arrays
275  integer,allocatable :: fs2full(:),indkk(:,:) !,fs2irr(:)
276  character(len=500) :: msg
277  type(fstab_t),pointer :: fs
278  type(t_tetrahedron) :: tetra
279 !arrays
280  !integer :: g0(3)
281  integer,allocatable :: full2ebands(:,:),bs2ibz(:)
282  real(dp) :: rlatt(3,3),klatt(3,3) !kibz(3),krot(3),
283  real(dp),allocatable :: kpt_full(:,:) !,fskpts(:,:)
284  real(dp), allocatable :: tmp_eigen(:),bdelta(:,:),btheta(:,:)
285 
286 ! *************************************************************************
287 
288  ABI_UNUSED(comm)
289 
290  !@fstab_t
291  call cwtime(cpu,wall,gflops,"start")
292 
293  if (any(cryst%symrel(:,:,1) /= identity_3d) .and. any(abs(cryst%tnons(:,1)) > tol10) ) then
294   MSG_ERROR('The first symmetry is not the identity operator!')
295  end if
296 
297  nkibz = ebands%nkpt
298 
299  !call kpts_ibz_from_kptrlatt(cryst, kptrlatt, ebands%kptopt, nshiftk, shiftk, nkibz, kibz, wtk, nkbz, kbz, &
300  ! new_kptrlatt, new_shiftk)  ! Optional
301 
302  ! Call smpbz to get the full grid of k-points `kpt_full`
303  ! brav1=1 is able to treat all bravais lattices (same option used in getkgrid)
304  mkpt= kptrlatt(1,1)*kptrlatt(2,2)*kptrlatt(3,3) &
305    +kptrlatt(1,2)*kptrlatt(2,3)*kptrlatt(3,1) &
306    +kptrlatt(1,3)*kptrlatt(2,1)*kptrlatt(3,2) &
307    -kptrlatt(1,2)*kptrlatt(2,1)*kptrlatt(3,3) &
308    -kptrlatt(1,3)*kptrlatt(2,2)*kptrlatt(3,1) &
309    -kptrlatt(1,1)*kptrlatt(2,3)*kptrlatt(3,2)
310 
311  ABI_STAT_MALLOC(kpt_full,(3,mkpt), ierr)
312  ABI_CHECK(ierr==0, 'allocating kpt_full')
313 
314  call smpbz(brav1,std_out,kptrlatt,mkpt,nkpt_full,nshiftk,option0,shiftk,kpt_full)
315 
316  ! Find correspondence BZ --> ebands%kpt
317  timrev = kpts_timrev_from_kptopt(ebands%kptopt)
318  sppoldbl = 1; if (any(cryst%symafm == -1) .and. ebands%nsppol==1) sppoldbl=2
319  ABI_MALLOC(indkk, (nkpt_full*sppoldbl,6))
320 
321  ! Compute k points from input file closest to the output file
322  call listkk(dksqmax,cryst%gmet,indkk,ebands%kptns,kpt_full,ebands%nkpt,nkpt_full,cryst%nsym,&
323     sppoldbl,cryst%symafm,cryst%symrel,timrev,use_symrec=.False.)
324 
325  if (dksqmax > tol12) then
326    write(msg, '(7a,es16.6,4a)' )&
327    'The WFK file cannot be used to start thee present calculation ',ch10,&
328    'It was asked that the wavefunctions be accurate, but',ch10,&
329    'at least one of the k points could not be generated from a symmetrical one.',ch10,&
330    'dksqmax=',dksqmax,ch10,&
331    'Action: check your WFK file and k-point input variables',ch10,&
332    '        (e.g. kptopt or shiftk might be wrong in the present dataset or the preparatory one.'
333    MSG_ERROR(msg)
334  end if
335 
336  call cwtime(cpu,wall,gflops,"stop")
337  write(msg,'(2(a,f8.2))')"fstab_init%listkk: cpu:",cpu,", wall: ",wall
338  call wrtout(std_out,msg,"COLL",do_flush=.True.)
339 
340  ABI_MALLOC(full2ebands, (6, nkpt_full))
341  full2ebands = 0
342 
343  do ik_full=1,nkpt_full
344    full2ebands(1, ik_full) = indkk(ik_full,1)     ! ik_ibz
345    full2ebands(2, ik_full) = indkk(ik_full,2)     ! isym
346    full2ebands(3, ik_full) = indkk(ik_full,6)     ! itimrev
347    full2ebands(4:6, ik_full) = indkk(ik_full,3:5) ! g0
348  end do
349 
350  call cwtime(cpu,wall,gflops,"start")
351 
352  ! Select only those k-points in the BZ close to the FS.
353  ABI_CHECK(fsewin > tol12, "fsewin < tol12")
354  elow = ebands%fermie - fsewin
355  ehigh = ebands%fermie + fsewin
356  ebis = elow - abs(elow) * 0.001_dp
357 
358  ! Allocate workspace arrays.
359  !ABI_MALLOC(fs2irr, (nkpt_full))
360  ABI_MALLOC(fs2full, (nkpt_full))
361 
362  do spin=1,ebands%nsppol
363    fs => fstab(spin)
364    ABI_MALLOC(fs%bstcnt_ibz, (2, nkibz))
365    fs%bstcnt_ibz = -1
366 
367    ! Find k-points on the FS(spin).
368    nkfs = 0
369    do ik_full=1,nkpt_full
370      ik_ibz = full2ebands(1, ik_full)
371      nband_k = ebands%nband(ik_ibz+(spin-1)*nkibz)
372 
373      blow = bisect(ebands%eig(:nband_k,ik_ibz,spin), ebis)
374      if (blow == 0) blow = 1
375      !if (blow == nband_k .or. blow == 0) cycle ! out of range
376      !write(std_out,*)"here with blow: ", blow,nband_k
377      !write(std_out,*)"eig_blow, eig_max, elow, ehigh:", &
378      !  ebands%eig(blow, ik_ibz, spin), ebands%eig(nband_k, ik_ibz, spin), elow,ehigh
379 
380      inwin = .False.; i1 = huge(1); i2 = -1
381      do band=blow,nband_k
382         !if (ebands%eig(band, ik_ibz, spin) > ehigh) exit
383         !write(std_out,*)band, ebands%eig(band, ik_ibz, spin) >= elow, ebands%eig(band, ik_ibz, spin) <= ehigh
384         if (ebands%eig(band, ik_ibz, spin) >= elow .and. ebands%eig(band, ik_ibz, spin) <= ehigh) then
385           inwin = .True.; i1 = min(i1, band); i2 = max(i2, band)
386         end if
387      end do
388 
389      if (inwin) then
390        ! Add this k-point and the corresponding bands.
391        !write(std_out,*)"in win"
392        nkfs = nkfs + 1
393        !fs2irr(nkfs) = ik_ibz
394        fs2full(nkfs) = ik_full
395        if (any(fs%bstcnt_ibz(:, ik_ibz) /= [-1, -1])) then
396          ABI_CHECK(all(fs%bstcnt_ibz(:, ik_ibz) == [i1, i2-i1+1]), "bstcnt_ibz!")
397        end if
398        fs%bstcnt_ibz(:, ik_ibz) = [i1, i2-i1+1]
399      end if
400    end do !ik_full
401 
402    ! @fstab_t
403    ! Build fstab_t for this spin.
404    fs%nkibz = nkibz; fs%nkfs = nkfs; fs%nktot = nkpt_full
405    ABI_MALLOC(fs%kpts, (3, nkfs))
406    ABI_MALLOC(fs%istg0, (6, nkfs))
407    do ik=1,nkfs
408      !ik_ibz = fs2irr(ik)
409      ik_full = fs2full(ik)
410      fs%kpts(:,ik) = kpt_full(:, ik_full)
411      fs%istg0(:, ik) = full2ebands(:, ik_full)
412    end do
413    fs%maxnb = maxval(fs%bstcnt_ibz(2, :))
414    call mkkptrank(fs%kpts,nkfs,fs%krank)
415  end do ! spin
416 
417  call cwtime(cpu,wall,gflops,"stop")
418  write(msg,'(2(a,f8.2))')"fstab_init%fs_build: cpu:",cpu,", wall: ",wall
419  call wrtout(std_out,msg,"COLL",do_flush=.True.)
420  call cwtime(cpu,wall,gflops,"start")
421 
422  ! fix window around fermie for tetrahedron or gaussian weight calculation
423  ! this is spin independent
424  nene = 100 ! TODO: make this variable and maybe temperature dependent???
425  deltaene = 2*fsewin/dble(nene-1)
426  ifermi = int(nene/2)
427  enemin = ebands%fermie - dble(ifermi-1)*deltaene
428  enemax = enemin + dble(nene-1)*deltaene
429 
430  ! Setup FS integration
431  do spin=1,ebands%nsppol
432    fs => fstab(spin)
433    fs%nene = nene
434    fs%enemin = enemin
435    fs%deltaene = deltaene
436    fs%nsig = 1
437    fs%integ_method = integ_method
438  end do
439 
440  if (integ_method == 2) then
441    rlatt = kptrlatt
442    call matr3inv(rlatt,klatt)
443 
444    ABI_MALLOC(bs2ibz, (nkpt_full))
445    bs2ibz = full2ebands(1, :)
446 
447    call init_tetra(bs2ibz, cryst%gprimd, klatt, kpt_full, nkpt_full, tetra, ierr, errstr)
448    ABI_CHECK(ierr==0, errstr)
449    ABI_FREE(bs2ibz)
450 
451    ABI_MALLOC(tmp_eigen, (nkibz))
452    ABI_MALLOC(btheta, (nene, nkibz))
453    ABI_MALLOC(bdelta, (nene,nkibz))
454 
455    max_occ = one
456    ! in spinor or spin polarized case, orbitals have occupation <= 1 instead of 2
457    !if (ebands%nsppol > 1) max_occ = one
458    ! this accounts for the doubling of the num of bands, even though spin channels are not well defined
459    !if (ebands%nspinor == 2) max_occ = half
460 
461    do spin=1,ebands%nsppol
462      fs => fstab(spin)
463      ABI_CALLOC(fs%tetra_wtk, (fs%maxnb, nkibz))
464      ABI_CALLOC(fs%tetra_wtk_ene, (fs%maxnb, nkibz, fs%nene))
465 
466      ! we have to pass full bands to get_tetra_weight
467      ! Here I create a full set of bands enclosing the states that will be used in the FS integration.
468      ! Then we have to rearrange the weights
469      bmin = huge(1); bmax = -huge(1)
470      do ik_ibz=1,nkibz
471        if (fs%bstcnt_ibz(1,ik_ibz) /= -1) then
472          bmin = min(bmin, fs%bstcnt_ibz(1,ik_ibz))
473        end if
474        if (fs%bstcnt_ibz(2,ik_ibz) /= -1) then
475          bmax = max(bmax, fs%bstcnt_ibz(1,ik_ibz) + fs%bstcnt_ibz(2,ik_ibz) - 1)
476        end if
477      end do
478      !write(std_out,*)"bmin, bmax for tetra: ",bmin,bmax
479      ABI_CHECK(bmin /= huge(1) .and. bmax /= -huge(1), "No point on the Fermi surface!")
480 
481      !call libtetrabz_dbldelta(ltetra, bvec, nb, nge, eig1, eig2, ngw, wght_bz, comm)
482 
483      do band=bmin,bmax
484        ! Get the contribution of this band
485        tmp_eigen = ebands%eig(band, :nkibz, spin)
486 
487        ! Calculate general integration weights at each irred kpoint 
488        ! as in Blochl et al PRB 49 16223 [[cite:Bloechl1994a]]
489        call tetra_blochl_weights(tetra,tmp_eigen,enemin,enemax,max_occ,fs%nene,nkibz,&
490          bcorr0,btheta,bdelta,xmpi_comm_self)
491 
492        do ik_ibz=1,nkibz
493          bstart_k = fs%bstcnt_ibz(1, ik_ibz); bstop_k = bstart_k + fs%bstcnt_ibz(2, ik_ibz) - 1
494          if (band >= bstart_k .and. band <= bstop_k) then
495            ! Save weights in the correct position.
496            ib = band - bstart_k + 1
497            fs%tetra_wtk(ib,ik_ibz) = bdelta(ifermi,ik_ibz) * nkibz
498            fs%tetra_wtk_ene(ib,ik_ibz,1:fs%nene) = bdelta(1:fs%nene,ik_ibz) * nkibz
499          end if
500        end do ! ik_ibz
501      end do ! band
502    end do ! sppol
503 
504    ABI_FREE(tmp_eigen)
505    ABI_FREE(btheta)
506    ABI_FREE(bdelta)
507 
508    call destroy_tetra(tetra)
509  end if
510 
511  !ABI_FREE(fs2irr)
512  ABI_FREE(fs2full)
513  ABI_FREE(kpt_full)
514  ABI_FREE(full2ebands)
515  ABI_FREE(indkk)
516 
517  call cwtime(cpu,wall,gflops,"stop")
518  write(msg,'(2(a,f8.2))')"fstab_init%fs_weights ",cpu,", wall: ",wall
519  call wrtout(std_out,msg,"COLL",do_flush=.True.)
520 
521 end subroutine fstab_init

m_fstab/fstab_print [ Functions ]

[ Top ] [ m_fstab ] [ Functions ]

NAME

  fstab_print

FUNCTION

  Print info on the object.

INPUTS

 [unit]=the unit number for output
 [prtvol]=verbosity level
 [mode_paral]=either "COLL" or "PERS"

OUTPUT

  Only printing.

PARENTS

      m_phgamma

CHILDREN

      wrtout

SOURCE

689 subroutine fstab_print(fstab, header, unit, prtvol, mode_paral)
690 
691 
692 !This section has been created automatically by the script Abilint (TD).
693 !Do not modify the following lines by hand.
694 #undef ABI_FUNC
695 #define ABI_FUNC 'fstab_print'
696 !End of the abilint section
697 
698  implicit none
699 
700 !Arguments ------------------------------------
701 !scalars
702  integer,optional,intent(in) :: prtvol,unit
703  character(len=4),optional,intent(in) :: mode_paral
704  character(len=*),optional,intent(in) :: header
705  type(fstab_t),target,intent(in) :: fstab(:)
706 
707 !Local variables-------------------------------
708 !scalars
709  integer :: my_unt,my_prtvol,spin
710  type(fstab_t),pointer :: fs
711  character(len=4) :: my_mode
712  character(len=500) :: msg
713 
714 ! *************************************************************************
715 
716  my_unt =std_out; if (present(unit)) my_unt = unit
717  my_prtvol=0    ; if (present(prtvol)) my_prtvol = prtvol
718  my_mode='COLL' ; if (present(mode_paral)) my_mode = mode_paral
719 
720  msg=' ==== Info on the fstab% object ==== '
721  if (PRESENT(header)) msg=' ==== '//TRIM(ADJUSTL(header))//' ==== '
722  call wrtout(my_unt,msg,my_mode)
723 
724  if (fstab(1)%integ_method == 1) then
725    write(std_out,"(a,i0)")"FS integration done with gaussians and nsig: ",fstab(1)%nsig
726  else if (fstab(1)%integ_method == 2) then
727    write(std_out,"(a)")"FS integration done with Tetrahedron method"
728  end if
729  write(std_out,"(a,i0)")"Total number of points in the full mesh: ",fstab(1)%nktot
730 
731  do spin=1,size(fstab)
732    fs => fstab(spin)
733    write(std_out,"(a,i0)")"For spin: ",spin
734    write(std_out,"(a,i0,a,f5.1,a)")&
735      "  Number of BZ k-points close to the Fermi surface: ",fs%nkfs," [",(100.0_dp*fs%nkfs)/fs%nktot," %]"
736    write(std_out,"(a,i0)")"  Maximum number of bands crossing the Fermi level: ",fs%maxnb
737    write(std_out,"(2(a,i0))")"  min band: ",minval(fs%bstcnt_ibz(1,:), mask=fs%bstcnt_ibz(1,:)/=-1)
738    write(std_out,"(2(a,i0))")"  Max band: ",maxval(fs%bstcnt_ibz(1,:)+fs%bstcnt_ibz(2,:)-1, mask=fs%bstcnt_ibz(1,:)/=-1)
739  end do
740 
741 end subroutine fstab_print

m_fstab/fstab_t [ Types ]

[ Top ] [ m_fstab ] [ Types ]

NAME

 fstab_t

FUNCTION

  Tables with the correspondence between points of the Fermi surface (FS) and the k-points in the
  IBZ (k-points found in ebands_t). We use `nsppol` fstab_t objects to account for spin polarization.

SOURCE

 60  type,public :: fstab_t
 61 
 62    integer :: nkfs
 63     ! Number of k-points on the Fermi-surface (FS-BZ).
 64 
 65    integer :: nktot
 66     ! Total number of k-points in the initial mesh.
 67 
 68    integer :: nkibz
 69     ! Number of points in the IBZ
 70 
 71    integer :: maxnb
 72    ! Max number of bands on the FS.
 73    ! TODO: Maybe maxnbfs
 74 
 75    ! real(dp) :: fermie
 76    ! Fermi energy
 77 
 78    integer :: integ_method
 79    ! Integration method. 1 for gaussian, 2 for tetrahedra
 80 
 81    integer :: nsig
 82    ! Number of smearing values used for Gaussian integration
 83 
 84    integer :: nene
 85    ! Number of chemical potential values used for inelastic integration
 86 
 87    real(dp) :: enemin
 88    ! Minimal chemical potential value used for inelastic integration
 89 
 90    real(dp) :: deltaene
 91    ! Chemical potential increment for inelastic integration
 92 
 93    type(kptrank_type) :: krank
 94    ! rank/inverse_rank pair for the k-points on the FS (kpts).
 95 
 96    integer,allocatable :: istg0(:,:)
 97    ! istg0(6, nkfs)
 98    ! Tables giving the correspondence between a point in the FS-BZ and the IBZ.
 99    !
100    ! istg0(1,:)      Mapping FS-BZ --> k-points in the IBZ (taken from ebands_t)
101    ! istg0(2,:)      The index of the symmetry S such that kfs = tim_sign * S(k_ibz) + G0
102    ! istg0(3,:)      1 if time-reversal was used to generate the k-point, 0 otherwise
103    ! istg0(4:6,:)    The reduced components of G0.
104 
105    integer,allocatable :: bstcnt_ibz(:,:)
106     ! bstcnt_ibz(2, nkibz)
107     ! The indices of the bands within the energy window (depends on fsk)
108     ! Note that we use the k-point index in the IBZ.
109     ! bstcnt(1, :)
110     ! The index of the first band inside the energy window (start)
111     ! bstcnt(2, :)
112     ! Number of bands on the FS (count)
113 
114    real(dp),allocatable :: kpts(:,:)
115    ! kpts(3,nkfs)
116    ! Reduced coordinates of the k-points on the Fermi surface.
117 
118    real(dp),allocatable :: tetra_wtk(:,:)
119    ! tetra_wtk(maxnb, nkibz)
120    ! Weights for FS integration with tetrahedron method
121    ! Note that the weights are dimensioned with nkibz
122 
123    real(dp),allocatable :: tetra_wtk_ene(:,:,:)
124    ! tetra_wtk_ene(maxnb, nkibz, nene)
125    ! Weights for FS integration with tetrahedron method
126    ! for all chemical potentials
127    ! Note that the weights are dimensioned with nkibz
128 
129  end type fstab_t
130 
131  public :: fstab_init            ! Initialize the object.
132  public :: fstab_free            ! Free memory.
133  public :: fstab_findkg0         ! Find the index of the k-point on the FS
134  public :: fstab_weights_ibz     ! Compute weights for FS integration.
135  public :: fstab_print           ! Print the object

m_fstab/fstab_weights_ibz [ Functions ]

[ Top ] [ m_fstab ] [ Functions ]

NAME

  fstab_weights_ibz

FUNCTION

  Return the weights for the integration on the Fermi-surface
  NB: single spin version - should have an array (nsppol) of the fstab objects

INPUTS

  ebands<ebands_type>=GS band structure.
  ik_ibz=Index of the k-point in the IBZ
  spin=Spin index
  sigmas

OUTPUT

   wtk(fs%nsig,fs%maxnb)=Weights for FS integration.

PARENTS

      m_ddk,m_phgamma

CHILDREN

      wrtout

SOURCE

602 subroutine fstab_weights_ibz(fs, ebands, ik_ibz, spin, sigmas, wtk, iene)
603 
604 
605 !This section has been created automatically by the script Abilint (TD).
606 !Do not modify the following lines by hand.
607 #undef ABI_FUNC
608 #define ABI_FUNC 'fstab_weights_ibz'
609 !End of the abilint section
610 
611  implicit none
612 
613 !Arguments ------------------------------------
614 !scalars
615  integer,intent(in) :: ik_ibz,spin
616  integer,intent(in),optional :: iene
617  type(fstab_t),intent(in) :: fs
618  type(ebands_t),intent(in) :: ebands
619 !arrays
620  real(dp),intent(in) :: sigmas(:) !fs%nsig)
621  real(dp),intent(out) :: wtk(fs%nsig,fs%maxnb)
622 
623 !Local variables-------------------------------
624 !scalars
625  integer :: ib,bstart_k,nband_k,band,isig
626  real(dp) :: arg
627  real(dp) :: chempot
628 
629 ! *************************************************************************
630 
631  !ik_ibz = fs%istg0(1, ik_bz)
632  bstart_k = fs%bstcnt_ibz(1, ik_ibz); nband_k = fs%bstcnt_ibz(2, ik_ibz)
633  ABI_CHECK(nband_k >= 1 .and. nband_k <= fs%maxnb, "wrong nband_k")
634 
635 ! TODO: add iene looping for chemical potential in gaussian case too
636  select case (fs%integ_method)
637  case (1)
638    chempot = ebands%fermie
639    if (present(iene)) then
640      chempot = fs%enemin + (iene-1)*fs%deltaene
641    end if
642    do ib=1,nband_k
643      band = ib + bstart_k - 1
644      arg = ebands%eig(band,ik_ibz,spin) - chempot
645      do isig=1,fs%nsig
646        wtk(isig,ib) = dirac_delta(arg, sigmas(isig))
647      end do
648    end do
649 
650  case (2)
651    if (present(iene)) then
652      wtk(1,1:nband_k) = fs%tetra_wtk_ene(1:nband_k, ik_ibz, iene)
653    else
654      wtk(1,1:nband_k) = fs%tetra_wtk(1:nband_k, ik_ibz)
655    end if
656 
657  case default
658    MSG_ERROR(sjoin("Wrong integration method:", itoa(fs%integ_method)))
659  end select
660 
661 end subroutine fstab_weights_ibz

m_fstab/mkqptequiv [ Functions ]

[ Top ] [ m_fstab ] [ Functions ]

NAME

 mkqptequiv

FUNCTION

 This routine determines the equivalence between
   1) qpoints and fermi surface kpoints
   2) qpoints under symmetry operations

INPUTS

   Cryst<crystal_t>=Info on unit cell and symmetries.
   kpt_phon = fermi surface kpoints
   nkpt_phon = number of kpoints in the full FS set
   nqpt = number of qpoints
   qpt_full = qpoint coordinates

OUTPUT

   FSfullpqtofull = mapping of k + q onto k' for k and k' in full BZ
   qpttoqpt(itim,isym,iqpt) = qpoint index which transforms to iqpt under isym and with time reversal itim.

NOTES

   REMOVED 3/6/2008: much too large matrix, and not used at present
       FStoqpt = mapping of kpoint pairs (1 irreducible and 1 full) to qpoints

PARENTS

      elphon,get_tau_k

CHILDREN

      destroy_kptrank,get_rank_1kpt,mkkptrank,wrtout

SOURCE

778 subroutine mkqptequiv(FSfullpqtofull,Cryst,kpt_phon,nkpt_phon,nqpt,qpttoqpt,qpt_full,mqtofull)
779 
780  use m_kptrank
781 
782 !This section has been created automatically by the script Abilint (TD).
783 !Do not modify the following lines by hand.
784 #undef ABI_FUNC
785 #define ABI_FUNC 'mkqptequiv'
786 !End of the abilint section
787 
788  implicit none
789 
790 !Arguments ------------------------------------
791 !scalars
792  integer,intent(in) :: nkpt_phon,nqpt
793  type(crystal_t),intent(in) :: Cryst
794 !arrays
795  integer,intent(out) :: FSfullpqtofull(nkpt_phon,nqpt),qpttoqpt(2,Cryst%nsym,nqpt)
796  integer,intent(out),optional :: mqtofull(nqpt)
797  real(dp),intent(in) :: kpt_phon(3,nkpt_phon),qpt_full(3,nqpt)
798 
799 !Local variables-------------------------------
800 !scalars
801  integer :: ikpt_phon,iFSqpt,iqpt,isym,symrankkpt_phon
802  character(len=500) :: message
803  type(kptrank_type) :: kptrank_t
804 !arrays
805  real(dp) :: tmpkpt(3),gamma_kpt(3)
806 
807 ! *************************************************************************
808 
809  call wrtout(std_out,' mkqptequiv : making rankkpt_phon and invrankkpt_phon',"COLL")
810 
811  call mkkptrank (kpt_phon,nkpt_phon,kptrank_t)
812 
813  FSfullpqtofull = -999
814  gamma_kpt(:) = zero
815 
816  do ikpt_phon=1,nkpt_phon
817    do iqpt=1,nqpt
818 !    tmpkpt = jkpt = ikpt + qpt
819      tmpkpt(:) = kpt_phon(:,ikpt_phon) + qpt_full(:,iqpt)
820 
821 !    which kpt is it among the full FS kpts?
822      call get_rank_1kpt (tmpkpt,symrankkpt_phon,kptrank_t)
823 
824      FSfullpqtofull(ikpt_phon,iqpt) = kptrank_t%invrank(symrankkpt_phon)
825      if (FSfullpqtofull(ikpt_phon,iqpt) == -1) then
826        MSG_ERROR("looks like no kpoint equiv to k+q !!!")
827      end if
828 
829    end do
830  end do
831 
832  if (present(mqtofull)) then
833    do iqpt=1,nqpt
834      tmpkpt(:) = gamma_kpt(:) - qpt_full(:,iqpt)
835 
836 !    which kpt is it among the full FS kpts?
837      call get_rank_1kpt (tmpkpt,symrankkpt_phon,kptrank_t)
838 
839      mqtofull(iqpt) = kptrank_t%invrank(symrankkpt_phon)
840      if (mqtofull(iqpt) == -1) then
841        MSG_ERROR("looks like no kpoint equiv to -q !!!")
842      end if
843    end do
844  end if
845 
846  call destroy_kptrank (kptrank_t)
847 
848 !start over with q grid
849  call wrtout(std_out,' mkqptequiv : FSfullpqtofull made. Do qpttoqpt',"COLL")
850 
851  call mkkptrank (qpt_full,nqpt,kptrank_t)
852 
853  qpttoqpt(:,:,:) = -1
854  do iFSqpt=1,nqpt
855    do isym=1,Cryst%nsym
856      tmpkpt(:) =  Cryst%symrec(:,1,isym)*qpt_full(1,iFSqpt) &
857 &     + Cryst%symrec(:,2,isym)*qpt_full(2,iFSqpt) &
858 &     + Cryst%symrec(:,3,isym)*qpt_full(3,iFSqpt)
859 
860      call get_rank_1kpt (tmpkpt,symrankkpt_phon,kptrank_t)
861      if (kptrank_t%invrank(symrankkpt_phon) == -1) then
862        message = "looks like no kpoint equiv to q by symmetry without time reversal!!!"
863        MSG_ERROR(message)
864      end if
865      qpttoqpt(1,isym,kptrank_t%invrank(symrankkpt_phon)) = iFSqpt
866 
867      tmpkpt = -tmpkpt
868      call get_rank_1kpt (tmpkpt,symrankkpt_phon,kptrank_t)
869      if (kptrank_t%invrank(symrankkpt_phon) == -1) then
870        message = ' mkqptequiv : Error : looks like no kpoint equiv to q by symmetry with time reversal!!!'
871        MSG_ERROR(message)
872      end if
873      qpttoqpt(2,isym,kptrank_t%invrank(symrankkpt_phon)) = iFSqpt
874    end do
875  end do
876 
877  call destroy_kptrank (kptrank_t)
878 
879 end subroutine mkqptequiv