TABLE OF CONTENTS
- ABINIT/m_fstab
- m_fstab/fstab_findkg0
- m_fstab/fstab_free
- m_fstab/fstab_get_dbldelta_weights
- m_fstab/fstab_init
- m_fstab/fstab_print
- m_fstab/fstab_t
ABINIT/m_fstab [ Modules ]
NAME
m_fstab
FUNCTION
Tools for the management of a set of Fermi surface k-points
COPYRIGHT
Copyright (C) 2008-2024 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 real(dp) :: abc(3) 597 598 ! ************************************************************************* 599 600 bstart_k = fs%bstart_cnt_ibz(1, ik_ibz); nband_k = fs%bstart_cnt_ibz(2, ik_ibz) 601 ABI_CHECK(nband_k >= 1 .and. nband_k <= fs%maxnb, "Wrong nband_k") 602 603 bstart_kq = fs%bstart_cnt_ibz(1, ikq_ibz); nband_kq = fs%bstart_cnt_ibz(2, ikq_ibz) 604 ABI_CHECK(nband_kq >= 1 .and. nband_kq <= fs%maxnb, "Wrong nband_kq") 605 606 wtk = zero 607 if (fs%eph_intmeth == 1 .or. nesting /= 0) then 608 ! Gaussian method: constant or adaptive method from group velocities if eph_fsmear is negative. 609 sigma = fs%eph_fsmear 610 use_adaptive = fs%eph_fsmear < zero .or. abs(fs%eph_intmeth) == 2 611 do ib2=1,nband_k 612 band2 = ib2 + bstart_k - 1 613 if (use_adaptive) then 614 !ori sigma = max(maxval([(abs(dot_product(fs%vk(:, ib2), fs%kmesh_cartvec(:,ii))), ii=1,3)]), fs%min_smear) 615 !workaround works with both ifort and ifx on oneapi 2024 616 !replace the implicit loop by an explicit one 617 do ii=1,3 618 abc(ii) = abs(dot_product(fs%vk(:, ib2), fs%kmesh_cartvec(:,ii))) 619 end do 620 sigma = max(maxval(abc), fs%min_smear) 621 !write(std_out, *)"sigma:", sigma * Ha_eV 622 end if 623 g2 = gaussian(ebands%eig(band2, ik_ibz, spin) - ebands%fermie, sigma) 624 do ib1=1,nband_kq 625 band1 = ib1 + bstart_kq - 1 626 if (use_adaptive) then 627 !ori sigma = max(maxval([(abs(dot_product(fs%vkq(:, ib1), fs%kmesh_cartvec(:,ii))), ii=1,3)]), fs%min_smear) 628 !replace the implicit loop by an explicit one 629 !workaround works with both ifort and ifx on oneapi 2024 630 do ii=1,3 631 abc(ii) = abs(dot_product(fs%vkq(:, ib1), fs%kmesh_cartvec(:,ii))) 632 end do 633 sigma = max(maxval(abc), fs%min_smear) 634 end if 635 g1 = gaussian(ebands%eig(band1, ikq_ibz, spin) - ebands%fermie, sigma) 636 wtk(ib1, ib2) = (g1 * g2) / fs%nktot 637 end do 638 end do 639 640 else if (abs(fs%eph_intmeth) == 2) then 641 ! Tetrahedron method. Copy weights in the correct position. 642 do ib2=1,nband_k 643 band2 = ib2 + bstart_k - fs%bmin 644 do ib1=1,nband_kq 645 band1 = ib1 + bstart_kq - fs%bmin 646 ! libtetrabz_dbldelta seems to report weights in this order. 647 wtk(ib1, ib2) = fs%dbldelta_tetra_weights_kfs(band1, band2, ik_fs) 648 end do 649 end do 650 651 else 652 ABI_ERROR(sjoin("Wrong integration method:", itoa(fs%eph_intmeth))) 653 end if 654 655 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
676 subroutine fstab_print(fstab, header, unit, prtvol) 677 678 !Arguments ------------------------------------ 679 !scalars 680 integer,optional,intent(in) :: prtvol,unit 681 character(len=*),optional,intent(in) :: header 682 class(fstab_t),target,intent(in) :: fstab(:) 683 684 !Local variables------------------------------- 685 !scalars 686 integer :: my_unt,my_prtvol,spin 687 class(fstab_t),pointer :: fs 688 character(len=5000) :: msg 689 690 ! ************************************************************************* 691 692 my_unt = std_out; if (present(unit)) my_unt = unit 693 my_prtvol = 0; if (present(prtvol)) my_prtvol = prtvol 694 695 msg = ' ==== Fermi surface info ==== ' 696 if (PRESENT(header)) msg=' ==== '//TRIM(ADJUSTL(header))//' ==== ' 697 write(my_unt, "(a)")trim(msg) 698 699 if (fstab(1)%eph_intmeth == 1) then 700 if (fstab(1)%eph_fsmear > zero) then 701 write(my_unt,"(a,f5.1,a)")" FS integration done with gaussian method and broadening:", & 702 fstab(1)%eph_fsmear * Ha_eV, " (meV)" 703 else 704 write(my_unt,"(a)")" FS integration done with adaptive gaussian method" 705 end if 706 else if (fstab(1)%eph_intmeth == 2) then 707 write(my_unt,"(a)")" FS integration done with tetrahedron method" 708 else if (fstab(1)%eph_intmeth == -2) then 709 write(my_unt,"(a)")" FS integration done with optimized tetrahedron method" 710 else 711 ABI_ERROR(sjoin("Invalid value for eph_intmeth:", itoa(fstab(1)%eph_intmeth))) 712 end if 713 714 write(my_unt,"(a,i0)")" Total number of k-points in the full mesh: ",fstab(1)%nktot 715 !write(my_unt,"(a,f5.1)")" Energy window: ",fstab(1)%eph_fsewin * Ha_eV, " (eV) 716 717 do spin=1,size(fstab) 718 fs => fstab(spin) 719 write(my_unt,"(a,i0)")" For spin: ",spin 720 write(my_unt,"(a,i0,a,f5.1,a)") & 721 " Number of BZ k-points close to the Fermi surface: ",fs%nkfs," [", (100.0_dp * fs%nkfs) / fs%nktot, " %]" 722 write(my_unt,"(a,i0)")" Maximum number of bands crossing the Fermi level: ",fs%maxnb 723 write(my_unt,"(2(a,i0))")" min band: ", minval(fs%bstart_cnt_ibz(1,:), mask=fs%bstart_cnt_ibz(1,:) /= -1) 724 write(my_unt,"(2(a,i0))")" Max band: ", maxval(fs%bstart_cnt_ibz(1,:) + fs%bstart_cnt_ibz(2,:) - 1, & 725 mask=fs%bstart_cnt_ibz(1,:) /= -1) 726 end do 727 728 end subroutine fstab_print
m_fstab/fstab_t [ 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