TABLE OF CONTENTS


ABINIT/m_fstab [ Modules ]

[ Top ] [ Modules ]

NAME

  m_fstab

FUNCTION

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

COPYRIGHT

  Copyright (C) 2008-2022 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 .

SOURCE

16 #if defined HAVE_CONFIG_H
17 #include "config.h"
18 #endif
19 
20 #include "abi_common.h"
21 
22 module m_fstab
23 
24  use defs_basis
25  use m_abicore
26  use m_xmpi
27  use m_errors
28  use m_krank
29  use m_htetra
30  use m_ebands
31  use m_crystal
32  use m_dtset
33 
34  use m_time,           only : cwtime, cwtime_report
35  use m_fstrings,       only : itoa, sjoin, ktoa
36  use m_numeric_tools,  only : bisect
37  use m_symtk,          only : matr3inv
38  use defs_datatypes,   only : ebands_t
39  use m_special_funcs,  only : gaussian
40  use m_kpts,           only : kpts_timrev_from_kptopt, smpbz, kpts_map
41 
42  implicit none
43 
44  private

m_fstab/fstab_findkg0 [ Functions ]

[ Top ] [ m_fstab ] [ Functions ]

NAME

  fstab_findkg0

FUNCTION

  Return the index `ik_fs` 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(:, ik_fs) + g0

SOURCE

539 integer function fstab_findkg0(fstab, kpt, g0) result(ik_fs)
540 
541 !Arguments ------------------------------------
542 !scalars
543  class(fstab_t),intent(in) :: fstab
544 !arrays
545  integer,intent(out) :: g0(3)
546  real(dp),intent(in) :: kpt(3)
547 
548 ! *************************************************************************
549 
550  ik_fs = fstab%krank%get_index(kpt)
551  if (ik_fs /= -1) then
552    g0 = nint(kpt - fstab%kpts(:, ik_fs))
553  else
554    g0 = huge(1)
555  end if
556 
557 end function fstab_findkg0

m_fstab/fstab_free [ Functions ]

[ Top ] [ m_fstab ] [ Functions ]

NAME

  fstab_free

FUNCTION

  Free memory

INPUTS

OUTPUT

SOURCE

192 subroutine fstab_free(fstab)
193 
194 !Arguments ------------------------------------
195  class(fstab_t),intent(inout) :: fstab
196 
197 ! ************************************************************************
198 
199  ! integer
200  ABI_SFREE(fstab%indkk_fs)
201  ABI_SFREE(fstab%bstart_cnt_ibz)
202 
203  ! real
204  ABI_SFREE(fstab%kpts)
205  ABI_SFREE(fstab%vk)
206  ABI_SFREE(fstab%vkq)
207  ABI_SFREE(fstab%tetra_wtk)
208  ABI_SFREE(fstab%tetra_wtk_ene)
209  ABI_SFREE(fstab%dbldelta_tetra_weights_kfs)
210 
211  ! types
212  call fstab%krank%free()
213 
214 end subroutine fstab_free

m_fstab/fstab_get_dbldelta_weights [ Functions ]

[ Top ] [ m_fstab ] [ Functions ]

NAME

  fstab_get_dbldelta_weights

FUNCTION

  Return the weights for the integration of the double-delta on the Fermi-surface

INPUTS

  ebands<ebands_type>=GS band structure.
  ik_ibz=Index of the k-point in the IBZ
  spin=Spin index
  nesting=Used in tetra. Set to 1 if tetra Weights cannot be computed with Tetrahedron due to nesting.
   In this case, we fallback to adaptive gaussian.

OUTPUT

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

SOURCE

581 subroutine fstab_get_dbldelta_weights(fs, ebands, ik_fs, ik_ibz, ikq_ibz, spin, nesting, wtk)
582 
583 !Arguments ------------------------------------
584 !scalars
585  integer,intent(in) :: ik_fs, ik_ibz, ikq_ibz, spin, nesting
586  class(fstab_t),intent(in) :: fs
587  type(ebands_t),intent(in) :: ebands
588 !arrays
589  real(dp),intent(out) :: wtk(fs%maxnb, fs%maxnb)
590 
591 !Local variables-------------------------------
592 !scalars
593  integer :: bstart_k, nband_k, bstart_kq, nband_kq, ib1, band1, ib2, band2, ii
594  logical :: use_adaptive
595  real(dp) :: g1, g2, sigma
596 
597 ! *************************************************************************
598 
599  bstart_k = fs%bstart_cnt_ibz(1, ik_ibz); nband_k = fs%bstart_cnt_ibz(2, ik_ibz)
600  ABI_CHECK(nband_k >= 1 .and. nband_k <= fs%maxnb, "Wrong nband_k")
601 
602  bstart_kq = fs%bstart_cnt_ibz(1, ikq_ibz); nband_kq = fs%bstart_cnt_ibz(2, ikq_ibz)
603  ABI_CHECK(nband_kq >= 1 .and. nband_kq <= fs%maxnb, "Wrong nband_kq")
604 
605  wtk = zero
606  if (fs%eph_intmeth == 1 .or. nesting /= 0) then
607    ! Gaussian method: constant or adaptive method from group velocities if eph_fsmear is negative.
608    sigma = fs%eph_fsmear
609    use_adaptive = fs%eph_fsmear < zero .or. abs(fs%eph_intmeth) == 2
610    do ib2=1,nband_k
611      band2 = ib2 + bstart_k - 1
612      if (use_adaptive) then
613        sigma = max(maxval([(abs(dot_product(fs%vk(:, ib2), fs%kmesh_cartvec(:,ii))), ii=1,3)]), fs%min_smear)
614        !write(std_out, *)"sigma:", sigma * Ha_eV
615      end if
616      g2 = gaussian(ebands%eig(band2, ik_ibz, spin) - ebands%fermie, sigma)
617      do ib1=1,nband_kq
618        band1 = ib1 + bstart_kq - 1
619        if (use_adaptive) then
620          sigma = max(maxval([(abs(dot_product(fs%vkq(:, ib1), fs%kmesh_cartvec(:,ii))), ii=1,3)]), fs%min_smear)
621        end if
622        g1 = gaussian(ebands%eig(band1, ikq_ibz, spin) - ebands%fermie, sigma)
623        wtk(ib1, ib2) = (g1 * g2) / fs%nktot
624      end do
625    end do
626 
627  else if (abs(fs%eph_intmeth) == 2) then
628    ! Tetrahedron method. Copy weights in the correct position.
629    do ib2=1,nband_k
630      band2 = ib2 + bstart_k - fs%bmin
631      do ib1=1,nband_kq
632        band1 = ib1 + bstart_kq - fs%bmin
633        ! libtetrabz_dbldelta seems to report weights in this order.
634        wtk(ib1, ib2) = fs%dbldelta_tetra_weights_kfs(band1, band2, ik_fs)
635      end do
636    end do
637 
638  else
639    ABI_ERROR(sjoin("Wrong integration method:", itoa(fs%eph_intmeth)))
640  end if
641 
642 end subroutine fstab_get_dbldelta_weights

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.
  dtset:
    eph_fsewin=Energy window in Hartree. Only states in [efermi-fsewin, efermi+fsewin] are included.
    eph_intmeth=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.

SOURCE

247 subroutine fstab_init(fstab, ebands, cryst, dtset, comm)
248 
249 !Arguments ------------------------------------
250 !scalars
251  integer,intent(in) :: comm
252  type(ebands_t),intent(in) :: ebands
253  type(crystal_t),intent(in) :: cryst
254  type(dataset_type),intent(in) :: dtset
255 !arrays
256  type(fstab_t),target,intent(out) :: fstab(ebands%nsppol)
257 
258 !Local variables-------------------------------
259 !scalars
260  integer,parameter :: option0 = 0, brav1 = 1, bcorr0 = 0
261  integer :: nkfs,spin,band,nband_k,i1,i2,ib,blow,ik_bz,ik_ibz,nkibz,timrev
262  integer :: ik,mkpt,nkbz,ierr, nene,ifermi
263  real(dp),parameter :: max_occ1 = one
264  real(dp) :: elow,ehigh,ebis,enemin,enemax,deltaene,cpu,wall,gflops
265  logical :: inwin
266  character(len=80) :: errstr
267  character(len=5000) :: msg
268  type(fstab_t),pointer :: fs
269  type(htetra_t) :: tetra
270  type(krank_t) :: krank
271 !arrays
272  integer :: kptrlatt(3,3)
273  integer,allocatable :: full2ebands(:,:),bz2ibz(:), fs2bz(:),indkk(:,:) !,fs2ibz(:)
274  real(dp),allocatable :: kbz(:,:), tmp_eigen(:),bdelta(:,:),btheta(:,:)
275  real(dp) :: rlatt(3,3), klatt(3,3)
276 
277 ! *************************************************************************
278 
279  call cwtime(cpu, wall, gflops, "start")
280 
281  if (any(cryst%symrel(:,:,1) /= identity_3d) .and. any(abs(cryst%tnons(:,1)) > tol10) ) then
282    ABI_ERROR('The first symmetry is not the identity operator!')
283  end if
284 
285  nkibz = ebands%nkpt
286  kptrlatt = dtset%kptrlatt
287  !call kpts_ibz_from_kptrlatt(cryst, kptrlatt, ebands%kptopt, dtset%nshiftk, dtset%shiftk, &
288  ! nkibz, kibz, wtk, nkbz, kbz, &
289  ! new_kptrlatt, new_shiftk)  ! Optional
290 
291  ! Call smpbz to get the full grid of k-points `kbz`
292  ! brav1=1 is able to treat all bravais lattices (same option used in getkgrid)
293  mkpt= kptrlatt(1,1)*kptrlatt(2,2)*kptrlatt(3,3) &
294    +kptrlatt(1,2)*kptrlatt(2,3)*kptrlatt(3,1) &
295    +kptrlatt(1,3)*kptrlatt(2,1)*kptrlatt(3,2) &
296    -kptrlatt(1,2)*kptrlatt(2,1)*kptrlatt(3,3) &
297    -kptrlatt(1,3)*kptrlatt(2,2)*kptrlatt(3,1) &
298    -kptrlatt(1,1)*kptrlatt(2,3)*kptrlatt(3,2)
299 
300  ABI_MALLOC(kbz, (3, mkpt))
301 
302  call smpbz(brav1, std_out, kptrlatt, mkpt, nkbz, dtset%nshiftk, option0, dtset%shiftk, kbz)
303 
304  ! Find correspondence BZ --> IBZ
305  ! Note that we use symrel so these tables can be used to symmetrize wavefunctions.
306  timrev = kpts_timrev_from_kptopt(ebands%kptopt)
307  ABI_MALLOC(indkk, (6, nkbz))
308 
309  krank = krank_from_kptrlatt(ebands%nkpt, ebands%kptns, kptrlatt, compute_invrank=.False.)
310 
311  if (kpts_map("symrel", timrev, cryst, krank, nkbz, kbz, indkk) /= 0) then
312    write(msg, '(10a)' ) &
313    'The WFK file cannot be used to start the present calculation ',ch10, &
314    'It was asked that the wavefunctions be accurate, but',ch10, &
315    'at least one of the k points could not be generated from a symmetrical one.',ch10, &
316    'Action: check your WFK file and k-point input variables',ch10, &
317    '        (e.g. kptopt or shiftk might be wrong in the present dataset or the preparatory one.'
318    ABI_ERROR(msg)
319  end if
320 
321  call krank%free()
322 
323  call cwtime_report(" fstab_init%krank", cpu, wall, gflops)
324 
325  ABI_MALLOC(full2ebands, (6, nkbz))
326  full2ebands = 0
327 
328  do ik_bz=1,nkbz
329    full2ebands(1, ik_bz) = indkk(1, ik_bz)      ! ik_ibz
330    full2ebands(2, ik_bz) = indkk(2, ik_bz)      ! isym
331    full2ebands(3:5, ik_bz) = indkk(3:5, ik_bz)  ! g0
332    full2ebands(6, ik_bz) = indkk(6, ik_bz)      ! itimrev
333  end do
334  ABI_FREE(indkk)
335 
336  ! Select only the k-points in the BZ that are sufficiently close to the FS.
337  ! FIXME: Do not know why but lambda depends on eph_fsewin if gaussian
338  ABI_CHECK(dtset%eph_fsewin > tol12, "dtset%eph_fsewin < tol12")
339  elow = ebands%fermie - dtset%eph_fsewin
340  ehigh = ebands%fermie + dtset%eph_fsewin
341  ebis = elow - abs(elow) * 0.001_dp
342 
343  ! Allocate workspace arrays.
344  !ABI_MALLOC(fs2ibz, (nkbz))
345  ABI_MALLOC(fs2bz, (nkbz))
346 
347  do spin=1,ebands%nsppol
348    fs => fstab(spin)
349    fs%spin = spin
350    ABI_MALLOC(fs%bstart_cnt_ibz, (2, nkibz))
351    fs%bstart_cnt_ibz = -1
352 
353    ! Find k-points on the FS associated to this spin.
354    nkfs = 0
355    do ik_bz=1,nkbz
356      ik_ibz = full2ebands(1, ik_bz)
357      nband_k = ebands%nband(ik_ibz + (spin-1)*nkibz)
358 
359      blow = bisect(ebands%eig(:nband_k, ik_ibz, spin), ebis)
360      if (blow == 0) blow = 1
361      !if (blow == nband_k .or. blow == 0) cycle ! out of range
362      !write(std_out,*)"here with blow: ", blow,nband_k
363      !write(std_out,*)"eig_blow, eig_max, elow, ehigh:", &
364      !                ebands%eig(blow, ik_ibz, spin), ebands%eig(nband_k, ik_ibz, spin), elow,ehigh
365 
366      inwin = .False.; i1 = huge(1); i2 = -1
367      do band=blow,nband_k
368         !if (ebands%eig(band, ik_ibz, spin) > ehigh) exit
369         !write(std_out,*)band, ebands%eig(band, ik_ibz, spin) >= elow, ebands%eig(band, ik_ibz, spin) <= ehigh
370         if (ebands%eig(band, ik_ibz, spin) >= elow .and. ebands%eig(band, ik_ibz, spin) <= ehigh) then
371           inwin = .True.; i1 = min(i1, band); i2 = max(i2, band)
372         end if
373      end do
374 
375      if (inwin) then
376        ! Add this k-point and the corresponding bands.
377        !write(std_out,*)"in win"
378        nkfs = nkfs + 1
379        !fs2ibz(nkfs) = ik_ibz
380        fs2bz(nkfs) = ik_bz
381        if (any(fs%bstart_cnt_ibz(:, ik_ibz) /= [-1, -1])) then
382          ABI_CHECK(all(fs%bstart_cnt_ibz(:, ik_ibz) == [i1, i2-i1+1]), "bstart_cnt_ibz!")
383        end if
384        fs%bstart_cnt_ibz(:, ik_ibz) = [i1, i2-i1+1]
385      end if
386    end do ! ik_bz
387 
388    ! Build fstab_t for this spin.
389    fs%nkibz = nkibz; fs%nkfs = nkfs; fs%nktot = nkbz
390    ABI_MALLOC(fs%kpts, (3, nkfs))
391    ABI_MALLOC(fs%indkk_fs, (6, nkfs))
392    do ik=1,nkfs
393      !ik_ibz = fs2ibz(ik)
394      ik_bz = fs2bz(ik)
395      fs%kpts(:,ik) = kbz(:, ik_bz)
396      fs%indkk_fs(:, ik) = full2ebands(:, ik_bz)
397    end do
398 
399    ! Define band indices enclosing states on the FS.
400    ! Note that we need all k-points for a given band when computing weights with tetrahedron.
401    ! This means that we have to be careful when selecting the weight associated to a given pair
402    ! (band_kq, kq), (band_k, k).
403    ! Then we have to rearrange the weights
404    fs%bmin = huge(1); fs%bmax = -huge(1)
405    do ik_ibz=1,nkibz
406      if (fs%bstart_cnt_ibz(1, ik_ibz) /= -1) then
407        fs%bmin = min(fs%bmin, fs%bstart_cnt_ibz(1,ik_ibz))
408      end if
409      if (fs%bstart_cnt_ibz(2, ik_ibz) /= -1) then
410        fs%bmax = max(fs%bmax, fs%bstart_cnt_ibz(1,ik_ibz) + fs%bstart_cnt_ibz(2,ik_ibz) - 1)
411      end if
412    end do
413 
414    !write(std_out,*)"bmin, bmax for tetra: ",fs%bmin, fs%bmax
415    ABI_CHECK(fs%bmin /= huge(1) .and. fs%bmax /= -huge(1), "No point on the Fermi surface!")
416    fs%maxnb = fs%bmax - fs%bmin + 1
417 
418    ! DEBUG: use same number of bands for each k-point on the FS.
419    !do ik_ibz=1,nkibz
420    !  if (fs%bstart_cnt_ibz(1, ik_ibz) /= -1) then
421    !    fs%bstart_cnt_ibz(1, ik_ibz) = fs%bmin
422    !    fs%bstart_cnt_ibz(2, ik_ibz) = fs%maxnb
423    !  end if
424    !end do
425 
426    ABI_CALLOC(fs%vk, (3, fs%maxnb))
427    ABI_CALLOC(fs%vkq, (3, fs%maxnb))
428 
429    fs%krank = krank_new(nkfs, fs%kpts)
430  end do ! spin
431 
432  call cwtime_report(" fstab_init%fs_build:", cpu, wall, gflops)
433 
434  ! fix window around fermie for tetrahedron or gaussian weight calculation
435  ! this is spin independent
436  nene = 100 ! TODO: make this variable and maybe temperature dependent???
437  deltaene = 2 * dtset%eph_fsewin / dble(nene-1)
438  ifermi = int(nene / 2)
439  enemin = ebands%fermie - dble(ifermi-1)*deltaene
440  enemax = enemin + dble(nene-1)*deltaene
441 
442  rlatt = kptrlatt
443  call matr3inv(rlatt, klatt)
444 
445  ! Setup FS integration
446  do spin=1,ebands%nsppol
447    fs => fstab(spin)
448    fs%nene = nene
449    fs%enemin = enemin
450    fs%deltaene = deltaene
451    fs%eph_intmeth = dtset%eph_intmeth
452    fs%eph_fsmear = dtset%eph_fsmear
453 
454    fs%klatt = klatt
455    fs%kmesh_cartvec(:, 1) = cryst%gprimd(:,1)*klatt(1,1) + cryst%gprimd(:,2)*klatt(2,1) + cryst%gprimd(:,3)*klatt(3,1)
456    fs%kmesh_cartvec(:, 2) = cryst%gprimd(:,1)*klatt(1,2) + cryst%gprimd(:,2)*klatt(2,2) + cryst%gprimd(:,3)*klatt(3,2)
457    fs%kmesh_cartvec(:, 3) = cryst%gprimd(:,1)*klatt(1,3) + cryst%gprimd(:,2)*klatt(2,3) + cryst%gprimd(:,3)*klatt(3,3)
458    ! TODO: It seems that two_pi is not needed here!
459    !fs%kmesh_cartvec = two_pi * fs%kmesh_cartvec
460    !do i1=1,3
461    !  write(std_out, *)"klatt:", klatt(:, i1)
462    !  write(std_out, *)"gprimd:", cryst%gprimd(:, i1)
463    !  write(std_out, *)"cartvec:", fs%kmesh_cartvec(:, i1)
464    !end do
465  end do
466 
467  if (abs(dtset%eph_intmeth) == 2) then
468    ! TODO: compute weights on the fly to reduce memory? nene should be set to zero if not used!
469    ABI_MALLOC(bz2ibz, (nkbz))
470    bz2ibz = full2ebands(1, :)
471 
472    call htetra_init(tetra, bz2ibz, cryst%gprimd, klatt, kbz, nkbz, ebands%kptns, nkibz, ierr, errstr, comm)
473    ABI_CHECK(ierr == 0, errstr)
474    ABI_FREE(bz2ibz)
475 
476    ABI_MALLOC(tmp_eigen, (nkibz))
477    ABI_MALLOC(btheta, (nene, nkibz))
478    ABI_MALLOC(bdelta, (nene, nkibz))
479 
480    do spin=1,ebands%nsppol
481      fs => fstab(spin)
482 
483      ! Allocate tables used to store tetrahedron weights.
484      ABI_CALLOC(fs%dbldelta_tetra_weights_kfs, (fs%maxnb, fs%maxnb, fs%nkfs))
485      ABI_CALLOC(fs%tetra_wtk, (fs%maxnb, nkibz))
486      ABI_CALLOC(fs%tetra_wtk_ene, (fs%maxnb, nkibz, fs%nene))
487 
488      do band=fs%bmin,fs%bmax
489        ! Get the contribution of this band
490        tmp_eigen = ebands%eig(band, :nkibz, spin)
491 
492        ! Calculate general integration weights at each irred kpoint
493        ! as in Blochl et al PRB 49 16223 [[cite:Bloechl1994a]]
494        call tetra%blochl_weights(tmp_eigen, enemin, enemax, max_occ1, fs%nene, nkibz, &
495          bcorr0, btheta, bdelta, xmpi_comm_self)
496 
497        ! Save weights in the correct position.
498        ib = band - fs%bmin + 1
499        do ik_ibz=1,nkibz
500          fs%tetra_wtk(ib, ik_ibz) = bdelta(ifermi, ik_ibz) * nkibz
501          fs%tetra_wtk_ene(ib, ik_ibz, 1:fs%nene) = bdelta(1:fs%nene, ik_ibz) * nkibz
502        end do
503      end do ! band
504    end do ! spin
505 
506    ABI_FREE(tmp_eigen)
507    ABI_FREE(btheta)
508    ABI_FREE(bdelta)
509    call tetra%free()
510  end if
511 
512  !ABI_FREE(fs2ibz)
513  ABI_FREE(fs2bz)
514  ABI_FREE(kbz)
515  ABI_FREE(full2ebands)
516 
517  call cwtime_report(" fstab_init%fs_weights:", cpu, wall, gflops)
518 
519 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

OUTPUT

  Only printing.

SOURCE

663 subroutine fstab_print(fstab, header, unit, prtvol)
664 
665 !Arguments ------------------------------------
666 !scalars
667  integer,optional,intent(in) :: prtvol,unit
668  character(len=*),optional,intent(in) :: header
669  class(fstab_t),target,intent(in) :: fstab(:)
670 
671 !Local variables-------------------------------
672 !scalars
673  integer :: my_unt,my_prtvol,spin
674  class(fstab_t),pointer :: fs
675  character(len=5000) :: msg
676 
677 ! *************************************************************************
678 
679  my_unt = std_out; if (present(unit)) my_unt = unit
680  my_prtvol = 0; if (present(prtvol)) my_prtvol = prtvol
681 
682  msg = ' ==== Fermi surface info ==== '
683  if (PRESENT(header)) msg=' ==== '//TRIM(ADJUSTL(header))//' ==== '
684  write(my_unt, "(a)")trim(msg)
685 
686  if (fstab(1)%eph_intmeth == 1) then
687    if (fstab(1)%eph_fsmear > zero) then
688      write(my_unt,"(a,f5.1,a)")" FS integration done with gaussian method and broadening:", &
689        fstab(1)%eph_fsmear * Ha_eV, " (meV)"
690    else
691      write(my_unt,"(a)")" FS integration done with adaptive gaussian method"
692    end if
693  else if (fstab(1)%eph_intmeth == 2) then
694    write(my_unt,"(a)")" FS integration done with tetrahedron method"
695  else if (fstab(1)%eph_intmeth == -2) then
696    write(my_unt,"(a)")" FS integration done with optimized tetrahedron method"
697  else
698    ABI_ERROR(sjoin("Invalid value for eph_intmeth:", itoa(fstab(1)%eph_intmeth)))
699  end if
700 
701  write(my_unt,"(a,i0)")" Total number of k-points in the full mesh: ",fstab(1)%nktot
702  !write(my_unt,"(a,f5.1)")" Energy window: ",fstab(1)%eph_fsewin * Ha_eV, " (eV)
703 
704  do spin=1,size(fstab)
705    fs => fstab(spin)
706    write(my_unt,"(a,i0)")" For spin: ",spin
707    write(my_unt,"(a,i0,a,f5.1,a)") &
708      "    Number of BZ k-points close to the Fermi surface: ",fs%nkfs," [", (100.0_dp * fs%nkfs) / fs%nktot, " %]"
709    write(my_unt,"(a,i0)")"    Maximum number of bands crossing the Fermi level: ",fs%maxnb
710    write(my_unt,"(2(a,i0))")"    min band: ", minval(fs%bstart_cnt_ibz(1,:), mask=fs%bstart_cnt_ibz(1,:) /= -1)
711    write(my_unt,"(2(a,i0))")"    Max band: ", maxval(fs%bstart_cnt_ibz(1,:) + fs%bstart_cnt_ibz(2,:) - 1, &
712                                                      mask=fs%bstart_cnt_ibz(1,:) /= -1)
713  end do
714 
715 end subroutine fstab_print

m_fstab/fstab_t [ Types ]

[ Top ] [ m_fstab ] [ Types ]

NAME

 fstab_t

FUNCTION

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

SOURCE

 58  type,public :: fstab_t
 59 
 60    integer :: spin
 61     ! Spin index
 62 
 63    integer :: nkfs
 64     ! Number of k-points on the Fermi-surface in the full BZ.
 65 
 66    integer :: nktot
 67     ! Total number of k-points in the initial mesh. Used to compute integrals in the BZ.
 68 
 69    integer :: nkibz
 70     ! Number of points in the IBZ
 71 
 72    integer :: bmin, bmax
 73     ! Min and max band index included in the calculation.
 74     ! Note that these values are obtained by taking the union over the k-points in the FS
 75     ! For each k-point, we usually have a different number of states crossing eF given by bstart_cnt_ibz
 76 
 77    integer :: maxnb
 78    ! Max number of bands on the FS (bmax - bmin + 1)
 79 
 80    integer :: eph_intmeth
 81    ! Integration method.
 82    ! 1 for gaussian (incuding adaptive broadening)
 83    ! |2| for tetrahedra.
 84    !     2 for the optimized tetrahedron method.
 85    !    -2 for the linear tetrahedron method.
 86    !
 87 
 88    integer :: nene
 89    ! Number of chemical potential values used for inelastic integration
 90 
 91    real(dp) :: eph_fsmear
 92    ! Gaussian broadening. Negative value activates adaptive gaussian broadening.
 93    ! See https://journals.aps.org/prb/pdf/10.1103/PhysRevB.92.075405
 94 
 95    real(dp) :: min_smear = tol9
 96    ! Used for the adaptive gaussian broadening: use min_smear if the broadening computed from the group velocity
 97    ! is smaller than this value to avoid divergences in the gaussian
 98 
 99    real(dp) :: enemin
100    ! Minimal chemical potential value used for inelastic integration
101 
102    real(dp) :: deltaene
103    ! Chemical potential increment for inelastic integration
104 
105    type(krank_t) :: krank
106    ! rank/inverse_rank pair for the k-points on the FS (kpts).
107 
108    integer,allocatable :: indkk_fs(:,:)
109    ! (6, nkfs)
110    ! Table giving the correspondence between a point in the FS-BZ and the IBZ:
111    !
112    !   indkk_fs(1,:)      Mapping FS-BZ --> k-points in the IBZ (taken from ebands_t)
113    !   indkk_fs(2,:)      The index of the symmetry S such that kfs = tim_sign * S(k_ibz) + G0
114    !   indkk_fs(3:5,:)    The reduced components of G0.
115    !   indkk_fs(6,:)      1 if time-reversal was used to generate the k-point, 0 otherwise
116 
117    integer,allocatable :: bstart_cnt_ibz(:,:)
118     ! (2, nkibz)
119     ! The indices of the bands within the energy window (depends on fsk)
120     ! Note that we use the k-point index in the IBZ.
121     !
122     !   bstcnt(1, :) The index of the first band inside the energy window (start)
123     !   bstcnt(2, :) Number of bands on the FS (count)
124 
125    real(dp) :: klatt(3, 3)
126     ! Reciprocal of lattice vectors for full kpoint grid. Used by init_tetra
127 
128    real(dp) :: kmesh_cartvec(3,3)
129     ! vectors defining the k-mesh (stored as column vector in Cartesian coords.
130     ! Used to implement the adaptive gaussian broadening.
131 
132    real(dp),allocatable :: kpts(:,:)
133    ! (3, nkfs)
134    ! Reduced coordinates of the BZ k-points on the Fermi surface.
135 
136    real(dp),allocatable :: vk(:,:), vkq(:,:)
137    ! (3, mnb)
138    ! Velocities in cartesian coordinates. Used to implement the adaptive gaussian broadening
139    ! Values are filled by the caller (e.g. phgamma) inside the loop over k-points.
140 
141    real(dp),allocatable :: tetra_wtk(:,:)
142    ! (maxnb, nkibz)
143    ! Weights for FS integration with tetrahedron method
144    ! Note that the weights are dimensioned with nkibz
145    ! (1, :) corresponds to %bmin
146 
147    real(dp),allocatable :: tetra_wtk_ene(:,:,:)
148    ! (maxnb, nkibz, nene)
149    ! Weights for FS integration with tetrahedron method for all chemical potentials
150    ! Note that the weights are dimensioned with nkibz
151    ! (1, :) corresponds to %bmin
152 
153    real(dp),allocatable :: dbldelta_tetra_weights_kfs(:,:,:)
154    ! (maxnb, maxnb, nkfs)
155    ! (1, 1, :) corresponds to %bmin
156 
157  contains
158 
159    procedure :: free => fstab_free
160    ! Free memory.
161 
162    procedure :: findkg0 => fstab_findkg0
163    ! Find the index of the k-point on the FS
164 
165    procedure :: get_dbldelta_weights => fstab_get_dbldelta_weights
166    ! Compute weights for the integration of the double-delta.
167 
168  end type fstab_t
169 
170  public :: fstab_init    ! Initialize the object.
171  public :: fstab_print   ! Print the object