TABLE OF CONTENTS
- ABINIT/m_sigtk
- m_sigmaph/sigtk_kpts_in_erange
- m_sigtk/sigtk_kcalc_from_erange
- m_sigtk/sigtk_kcalc_from_gaps
- m_sigtk/sigtk_kcalc_from_nkptgw
- m_sigtk/sigtk_kcalc_from_qprange
- m_sigtk/sigtk_sigma_tables
ABINIT/m_sigtk [ Modules ]
NAME
m_sigtk
FUNCTION
Helper functions common to electron self-energy calculations. Provides tools to: - Define list of k-points and bands in sel-energy matrix elements from input variables.
COPYRIGHT
Copyright (C) 2008-2024 ABINIT group (MG) 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
18 #if defined HAVE_CONFIG_H 19 #include "config.h" 20 #endif 21 22 #include "abi_common.h" 23 24 module m_sigtk 25 26 use defs_basis 27 use m_abicore 28 use m_errors 29 use m_ebands 30 use m_crystal 31 use m_xmpi 32 #ifdef HAVE_NETCDF 33 use netcdf 34 #endif 35 use m_nctk 36 use m_hdr 37 use m_dtset 38 use m_krank 39 40 use m_build_info, only : abinit_version 41 use m_fstrings, only : sjoin, ltoa, strcat, itoa, ftoa 42 use m_io_tools, only : open_file 43 use defs_datatypes, only : ebands_t, pseudopotential_type 44 use defs_wvltypes, only : wvl_internal_type 45 use m_pawtab, only : pawtab_type 46 use m_kpts, only : kpts_ibz_from_kptrlatt, kpts_timrev_from_kptopt, kpts_map 47 48 implicit none 49 50 private 51 52 public :: sigtk_kcalc_from_nkptgw 53 public :: sigtk_kcalc_from_qprange 54 public :: sigtk_kcalc_from_gaps 55 public :: sigtk_kcalc_from_erange 56 public :: sigtk_kpts_in_erange 57 public :: sigtk_sigma_tables
m_sigmaph/sigtk_kpts_in_erange [ Functions ]
[ Top ] [ m_sigmaph ] [ Functions ]
NAME
sigtk_kpts_in_erange
FUNCTION
Use star functions interpolation and [[einterp]] to interpolate KS energies onto dense k-mesh defined by [[sigma_ngkpt]] and [[sigma_shiftk]]. find k-points inside (electron/hole) pockets according to the values specifed by [[sigma_erange]]. write kerange.nc file with the tables required by abinit to automate nscf band structure calculations mainly used to prepare eph calculations in which only selected k-points are nededed (imaginary part of self-energies).
INPUTS
dtset <dataset_type>=all input variables for this dataset cryst<crystal_t>=Crystalline structure ebands<ebands_t>=The GS KS band structure (energies, occupancies, k-weights...) psps <pseudopotential_type>=all the information about psps pawtab(ntypat*usepaw) <type(pawtab_type)>=paw tabulated starting data prefix=Prefix for output file. comm: MPI communicator.
SOURCE
559 subroutine sigtk_kpts_in_erange(dtset, cryst, ebands, psps, pawtab, prefix, comm) 560 561 !Arguments ------------------------------------ 562 !scalars 563 type(dataset_type),intent(in) :: dtset 564 type(crystal_t),intent(in) :: cryst 565 type(ebands_t),intent(in) :: ebands 566 type(pseudopotential_type),intent(in) :: psps 567 character(len=*),intent(in) :: prefix 568 integer,intent(in) :: comm 569 !arrays 570 type(pawtab_type),intent(in) :: pawtab(dtset%ntypat*psps%usepaw) 571 572 !Local variables ------------------------------ 573 !scalars 574 integer,parameter :: master = 0, pertcase0 = 0, image1 = 1 575 integer :: ii, my_rank, nprocs, spin, ikf_ibz, band, nkpt_inerange, gap_err, unt, ncid, cnt, ncerr 576 logical :: assume_gap 577 real(dp) :: ee, cmin, vmax 578 character(len=500) :: msg 579 character(len=fnlen) :: path 580 type(ebands_t) :: fine_ebands 581 type(gaps_t) :: gaps, fine_gaps 582 type(wvl_internal_type) :: dummy_wvl 583 type(hdr_type) :: fine_hdr 584 !arrays 585 integer :: fine_kptrlatt(3,3), band_block(2), unts(2) 586 integer,allocatable :: kshe_mask(:,:,:), krange2ibz(:) 587 real(dp) :: params(4) 588 589 ! ************************************************************************* 590 591 my_rank = xmpi_comm_rank(comm); nprocs = xmpi_comm_size(comm) 592 593 ! (-num, -num) activate treatment of metals with energy window around Efermi. 594 assume_gap = .not. all(dtset%sigma_erange < zero) 595 unts = [std_out, ab_out] 596 597 if (my_rank == master) then 598 call wrtout(unts, sjoin(ch10, repeat("=", 92))) 599 call wrtout(unts, " Using SKW interpolation to interpolate KS energies onto dense k-mesh.") 600 call wrtout(unts, sjoin(" defined by sigma_ngkpt:", trim(ltoa(dtset%sigma_ngkpt)))) 601 ABI_CHECK(allocated(dtset%sigma_shiftk), "sigma_nshiftk must be specified in input.") 602 write(std_out, "(2a)") " and sigma_shiftk shifts:" 603 do ii=1,dtset%nshiftk 604 call wrtout(unts, sjoin(itoa(ii), ltoa(dtset%sigma_shiftk(:, ii)))) 605 end do 606 607 if (assume_gap) then 608 call wrtout(unts, " Finding k-points inside (electron/hole) pockets (assuming semiconductor).") 609 else 610 call wrtout(unts, " Finding k-points inside energy window around Fermi level (assuming metal).") 611 end if 612 write(msg, "(a, 2(f6.3, 1x), a)")" Using sigma_erange: ", dtset%sigma_erange(:) * Ha_eV, " (eV)" 613 call wrtout(unts, msg) 614 call wrtout(unts, sjoin(" SKW parameters (einterp): ", ltoa(dtset%einterp))) 615 call wrtout(unts, sjoin(repeat("=", 92), ch10)) 616 !call ebands_print(ebands, header, unit=std_out, prtvol=dtset%prtvol) 617 618 ! Consistency check. 619 if (all(dtset%sigma_erange == zero)) then 620 ABI_ERROR("sigma_erange must be specified in input when calling sigtk_kpts_in_erange.") 621 end if 622 if (all(dtset%sigma_ngkpt == 0)) then 623 ABI_ERROR("sigma_ngkpt must be specified in input when calling sigtk_kpts_in_erange.") 624 end if 625 end if 626 627 if (assume_gap) then 628 ! Compute gaps using input ebands. 629 gaps = ebands_get_gaps(ebands, gap_err) 630 if (gap_err /= 0) then 631 ABI_ERROR("Cannot compute fundamental and direct gap (likely metal).") 632 end if 633 634 if (my_rank == master) then 635 call gaps%print(header="Gaps from input WFK", unit=std_out) 636 call gaps%print(header="Gaps from input WFK", unit=ab_out) 637 end if 638 call gaps%free() 639 else 640 call wrtout(unts, sjoin("Using Fermi level:", ftoa(ebands%fermie * Ha_eV, fmt="f6.2"), " (eV)")) 641 end if 642 643 ! Interpolate band energies with star functions. 644 ! In the EPH code, we will need eigens in the IBZ to compute efermi not just energies inside pockets. 645 fine_kptrlatt = 0 646 do ii=1,3 647 fine_kptrlatt(ii, ii) = dtset%sigma_ngkpt(ii) 648 end do 649 band_block = [1, ebands%mband] 650 params = 0; params(1) = 1; params(2) = 5; if (nint(dtset%einterp(1)) == 1) params = dtset%einterp 651 652 fine_ebands = ebands_interp_kmesh(ebands, cryst, params, fine_kptrlatt, & 653 dtset%sigma_nshiftk, dtset%sigma_shiftk, band_block, comm) 654 fine_ebands%istwfk = 1 655 656 call ebands_update_occ(fine_ebands, dtset%spinmagntarget, prtvol=dtset%prtvol) 657 call ebands_print(fine_ebands, header="FINE EBANDS", unit=std_out, prtvol=dtset%prtvol) 658 659 if (assume_gap) then 660 ! Compute gaps using fine_ebands. 661 fine_gaps = ebands_get_gaps(fine_ebands, gap_err) 662 if (gap_err /= 0) then 663 ABI_ERROR("Cannot compute fundamental and direct gap (likely metal).") 664 end if 665 666 if (my_rank == master) then 667 call fine_gaps%print(header="Gaps from SKW interpolated eigenvalues", unit=std_out) 668 call fine_gaps%print(header="Gaps from SKW interpolated eigenvalues", unit=ab_out) 669 end if 670 end if 671 672 ! Build new header with fine k-mesh (note kptrlatt_orig == kptrlatt) 673 call hdr_init_lowlvl(fine_hdr, fine_ebands, psps, pawtab, dummy_wvl, abinit_version, pertcase0, & 674 dtset%natom, dtset%nsym, dtset%nspden, dtset%ecut, dtset%pawecutdg, dtset%ecutsm, dtset%dilatmx, & 675 dtset%intxc, dtset%ixc, dtset%stmbias, dtset%usewvl, dtset%pawcpxocc, dtset%pawspnorb, dtset%ngfft, dtset%ngfftdg, & 676 dtset%so_psp, dtset%qptn, cryst%rprimd, cryst%xred, cryst%symrel, cryst%tnons, cryst%symafm, cryst%typat, & 677 dtset%amu_orig(:, image1), dtset%icoulomb, & 678 dtset%kptopt, dtset%nelect, dtset%ne_qFD, dtset%nh_qFD, dtset%ivalence, dtset%cellcharge(1), & 679 fine_kptrlatt, fine_kptrlatt, dtset%sigma_nshiftk, dtset%sigma_nshiftk, dtset%sigma_shiftk, dtset%sigma_shiftk) 680 681 ! Find k-points inside sigma_erange energy window. 682 ! Set entry to the number of states inside the pocket at (ikpt, spin) 683 ! (last index discerns between hole and electron pockets) 684 ABI_ICALLOC(kshe_mask, (fine_ebands%nkpt, ebands%nsppol, 2)) 685 686 do spin=1,ebands%nsppol 687 ! Get CBM and VBM with some tolerance. 688 if (assume_gap) then 689 vmax = fine_gaps%vb_max(spin) + tol2 * eV_Ha 690 cmin = fine_gaps%cb_min(spin) - tol2 * eV_Ha 691 else 692 ! Note that we use the Fermi level from ebands instead of fine_ebands. 693 vmax = ebands%fermie 694 cmin = ebands%fermie 695 end if 696 697 do ikf_ibz=1,fine_ebands%nkpt 698 do band=1,ebands%mband 699 ee = fine_ebands%eig(band, ikf_ibz, spin) 700 ! Check whether the interpolated eigenvalue is inside the sigma_erange window. 701 if (abs(dtset%sigma_erange(1)) > zero) then 702 if (ee <= vmax .and. vmax - ee <= abs(dtset%sigma_erange(1))) then 703 kshe_mask(ikf_ibz, spin, 1) = kshe_mask(ikf_ibz, spin, 1) + 1; exit 704 end if 705 end if 706 if (abs(dtset%sigma_erange(2)) > zero) then 707 if (ee >= cmin .and. ee - cmin <= abs(dtset%sigma_erange(2))) then 708 kshe_mask(ikf_ibz, spin, 2) = kshe_mask(ikf_ibz, spin, 2) + 1; exit 709 end if 710 end if 711 end do 712 end do 713 end do 714 715 ! Build list of k-points inside pockets. Use over dimensioned array. 716 cnt = count(kshe_mask /= 0) 717 ABI_MALLOC(krange2ibz, (cnt)) 718 cnt = 0 719 do ikf_ibz=1,fine_ebands%nkpt 720 if (any(kshe_mask(ikf_ibz,:,:) /= 0)) then 721 cnt = cnt + 1 722 krange2ibz(cnt) = ikf_ibz 723 end if 724 end do 725 nkpt_inerange = cnt 726 727 ! Possible extensions that may be implemented at this level: 728 ! 1. Find image points in the BZ? 729 ! 2. Compute tetra and q-points for EPH calculation or use +/- wmax window and heuristic approach in sigmaph at runtime? 730 ! 3. Compute SKW 1st and 2nd derivatives needed to treat Frohlich? 731 732 ! Write output files with k-point list. 733 if (my_rank == master .and. len_trim(prefix) /= 0) then 734 write(std_out, "(a,i0,a,f5.1,a)")" Found: ", nkpt_inerange, " kpoints in sigma_erange energy windows. (nkeff / nkibz): ", & 735 (100.0_dp * nkpt_inerange) / fine_ebands%nkpt, " [%]" 736 737 ! Write text file with Abinit input variables (mainly for testing purposes). 738 path = strcat(prefix, "_KERANGE") 739 if (open_file(path, msg, newunit=unt, form="formatted") /= 0) then 740 ABI_ERROR(msg) 741 end if 742 write(unt, "(a)")"kptopt 0" 743 write(unt, "(a, i0)")"nkpt ", nkpt_inerange 744 write(unt, "(a)")"kpt" 745 do ii=1,nkpt_inerange 746 write(unt, "(3(es16.8,1x))") fine_ebands%kptns(:, krange2ibz(ii)) 747 end do 748 write(unt, "(a, i0)")"wtk" 749 do ii=1,nkpt_inerange 750 write(unt, "(es16.8)") fine_ebands%wtk(krange2ibz(ii)) 751 end do 752 close(unt) 753 754 ! Write netcdf file used to perform NSCF run and EPH calculations with eph_task = -4. 755 path = strcat(prefix, "_KERANGE.nc") 756 #ifdef HAVE_NETCDF 757 NCF_CHECK(nctk_open_create(ncid, path, xmpi_comm_self)) 758 ! Write crystalline structure, fine_hdr and fine_ebands defined on the fine k-mesh. 759 ! fine_ebands will be used to compare with the ab-initio NSCF eigenvalues. 760 ! 761 ! TODO: The size of the KERANGE.nc quickly increases with the k-mesh. 762 ! It is ~700 Mb for a ~ 300^3 grid due to occ and eigens 763 ! But these quantities are now used in inkpts so it may be possible to avoid writing them to disk. 764 ! 765 NCF_CHECK(fine_hdr%ncwrite(ncid, fform_from_ext("KERANGE.nc"), nc_define=.True.)) 766 NCF_CHECK(cryst%ncwrite(ncid)) 767 NCF_CHECK(ebands_ncwrite(fine_ebands, ncid)) 768 NCF_CHECK(nctk_def_dims(ncid, [nctkdim_t("nkpt_inerange", nkpt_inerange)], defmode=.True.)) 769 ! Define extra arrays. 770 ncerr = nctk_def_arrays(ncid, [ & 771 nctkarr_t("kshe_mask", "int", "number_of_kpoints, number_of_spins, two"), & 772 nctkarr_t("krange2ibz", "int", "nkpt_inerange"), & 773 nctkarr_t("sigma_erange", "dp", "two"), & 774 nctkarr_t("einterp", "dp", "four") & 775 ], defmode=.True.) 776 NCF_CHECK(ncerr) 777 ! Write extra arrays. 778 NCF_CHECK(nctk_set_datamode(ncid)) 779 NCF_CHECK(nf90_put_var(ncid, nctk_idname(ncid, "kshe_mask"), kshe_mask)) 780 NCF_CHECK(nf90_put_var(ncid, nctk_idname(ncid, "krange2ibz"), krange2ibz(1:nkpt_inerange))) 781 NCF_CHECK(nf90_put_var(ncid, nctk_idname(ncid, "sigma_erange"), dtset%sigma_erange)) 782 NCF_CHECK(nf90_put_var(ncid, nctk_idname(ncid, "einterp"), params)) 783 NCF_CHECK(nf90_close(ncid)) 784 #endif 785 end if 786 787 ABI_FREE(kshe_mask) 788 ABI_FREE(krange2ibz) 789 790 call fine_gaps%free() 791 call ebands_free(fine_ebands) 792 call fine_hdr%free() 793 794 end subroutine sigtk_kpts_in_erange
m_sigtk/sigtk_kcalc_from_erange [ Functions ]
[ Top ] [ m_sigtk ] [ Functions ]
NAME
sigtk_kcalc_from_erange
FUNCTION
Select list of k-points and bands for self-energy matrix elements on the basis of their positions wrt to the (band edges|fermi level) and the value of sigma_erange. Useful when computing electron-lifetimes for transport calculations. INPUT dtset<dataset_type>=All input variables for this dataset. cryst<crystal_t>=Crystalline structure ebands<ebands_t>=The GS KS band structure (energies, occupancies, k-weights...) gaps<gaps_t>=Store location of gaps. comm: MPI communicator.
OUTPUT
nkcalc: Number of k-points in self-energy matrix elements. kcalc(3, nkcalc): List of k-points where the self-energy is computed. bstart_ks(nkcalc, nsppol): Initial KS band index included in self-energy matrix elements for each k-point in kcalc. nbcalc_ks(nkcalc, nsppol): Number of bands included in self-energy matrix elements for each k-point in kcalc.
SOURCE
354 subroutine sigtk_kcalc_from_erange(dtset, cryst, ebands, gaps, nkcalc, kcalc, bstart_ks, nbcalc_ks, comm) 355 356 !Arguments ------------------------------------ 357 type(dataset_type),intent(in) :: dtset 358 type(crystal_t),intent(in) :: cryst 359 type(ebands_t),intent(in) :: ebands 360 type(gaps_t),intent(in) :: gaps 361 integer,intent(in) :: comm 362 integer,intent(out) :: nkcalc 363 !arrays 364 real(dp),allocatable,intent(out) :: kcalc(:,:) 365 integer,allocatable,intent(out) :: bstart_ks(:,:) 366 integer,allocatable,intent(out) :: nbcalc_ks(:,:) 367 368 !Local variables ------------------------------ 369 !scalars 370 integer,parameter :: master = 0 371 integer :: spin, ik, band, ii, ic, nsppol, tmp_nkpt, timrev, sigma_nkbz, my_rank 372 logical :: found 373 real(dp) :: cmin, vmax, ee 374 logical :: assume_gap 375 character(len=500) :: msg 376 type(krank_t) :: krank 377 !arrays 378 integer :: kptrlatt(3,3), unts(1) 379 integer,allocatable :: ib_work(:,:,:), sigmak2ebands(:), indkk(:,:) 380 integer :: kpos(ebands%nkpt) 381 real(dp),allocatable :: sigma_wtk(:),sigma_kbz(:,:),tmp_kcalc(:,:) 382 383 ! ************************************************************************* 384 385 my_rank = xmpi_comm_rank(comm) !; nprocs = xmpi_comm_size(comm) 386 unts = [std_out] 387 assume_gap = .not. all(dtset%sigma_erange < zero) 388 389 if (my_rank == master) then 390 write(std_out, "(a)")" Selecting k-points and bands according to their position wrt the band edges (sigma_erange)." 391 write(std_out, "(a, 2(f6.3, 1x), a)")" sigma_erange: ", dtset%sigma_erange(:) * Ha_eV, " (eV)" 392 if (assume_gap) then 393 call gaps%print(unit=std_out) 394 ABI_CHECK(maxval(gaps%ierr) == 0, "sigma_erange 0 cannot be used because I cannot find the gap (gap_err !=0)") 395 end if 396 end if 397 398 if (any(dtset%sigma_ngkpt /= 0)) then 399 call wrtout(std_out, sjoin(" Generating initial list of k-points from sigma_nkpt:", ltoa(dtset%sigma_ngkpt))) 400 ! Get tentative tmp_nkpt and tmp_kcalc from sigma_ngkpt. 401 kptrlatt = 0 402 kptrlatt(1,1) = dtset%sigma_ngkpt(1); kptrlatt(2,2) = dtset%sigma_ngkpt(2); kptrlatt(3,3) = dtset%sigma_ngkpt(3) 403 call kpts_ibz_from_kptrlatt(cryst, kptrlatt, dtset%kptopt, dtset%sigma_nshiftk, dtset%sigma_shiftk, & 404 tmp_nkpt, tmp_kcalc, sigma_wtk, sigma_nkbz, sigma_kbz) 405 406 ABI_FREE(sigma_kbz) 407 ABI_FREE(sigma_wtk) 408 409 ! Map tmp_kcalc to ebands%kpts 410 timrev = kpts_timrev_from_kptopt(ebands%kptopt) 411 412 ABI_MALLOC(indkk, (6, tmp_nkpt)) 413 414 krank = krank_from_kptrlatt(ebands%nkpt, ebands%kptns, ebands%kptrlatt, compute_invrank=.False.) 415 416 if (kpts_map("symrec", timrev, cryst, krank, tmp_nkpt, tmp_kcalc, indkk) /= 0) then 417 write(msg, '(3a)' )& 418 "At least one of the k-points could not be generated from a symmetrical one in the WFK.",ch10,& 419 'Action: check your WFK file and the value of sigma_nkpt, sigma_shiftk in the input file.' 420 ABI_ERROR(msg) 421 end if 422 423 call krank%free() 424 425 ABI_MALLOC(sigmak2ebands, (tmp_nkpt)) 426 sigmak2ebands = indkk(1, :) 427 ABI_FREE(tmp_kcalc) 428 ABI_FREE(indkk) 429 430 else 431 ! Include all the k-points in the IBZ in the initial list. 432 call wrtout(std_out, " Generating initial list of k-points from input ebands%kptns.") 433 tmp_nkpt = ebands%nkpt 434 ! Trivial map 435 ABI_MALLOC(sigmak2ebands, (tmp_nkpt)) 436 sigmak2ebands = [(ii, ii=1, ebands%nkpt)] 437 end if 438 439 nsppol = ebands%nsppol 440 ABI_MALLOC(ib_work, (2, tmp_nkpt, nsppol)) 441 442 do spin=1,nsppol 443 444 if (assume_gap) then 445 ! Get CBM and VBM with some tolerance 446 vmax = gaps%vb_max(spin) + tol2 * eV_Ha 447 cmin = gaps%cb_min(spin) - tol2 * eV_Ha 448 else 449 vmax = ebands%fermie 450 cmin = ebands%fermie 451 end if 452 453 do ii=1,tmp_nkpt 454 ! Index of k-point in ebands. 455 ik = sigmak2ebands(ii) 456 ! Will use this initial values to understand if k-point is in energy window. 457 ib_work(1, ii, spin) = huge(1) 458 ib_work(2, ii, spin) = -huge(1) 459 do band=1,ebands%nband(ik + (spin-1) * ebands%nkpt) 460 ee = ebands%eig(band, ik, spin) 461 if (abs(dtset%sigma_erange(1)) > zero) then 462 if (ee <= vmax .and. vmax - ee <= abs(dtset%sigma_erange(1))) then 463 ib_work(1, ii, spin) = min(ib_work(1, ii, spin), band) 464 ib_work(2, ii, spin) = max(ib_work(2, ii, spin), band) 465 !write(std_out, *), "Adding valence band", band, " with ee [eV]: ", ee * Ha_eV 466 end if 467 end if 468 if (abs(dtset%sigma_erange(2)) > zero) then 469 if (ee >= cmin .and. ee - cmin <= abs(dtset%sigma_erange(2))) then 470 ib_work(1, ii, spin) = min(ib_work(1, ii, spin), band) 471 ib_work(2, ii, spin) = max(ib_work(2, ii, spin), band) 472 !write(std_out, *)"Adding conduction band", band, " with ee [eV]: ", ee * Ha_eV 473 end if 474 end if 475 end do 476 end do 477 end do 478 479 ! Now we can define the list of k-points and the bands range. 480 ! The main problem here is that kptgw and nkptgw do not depend on the spin and therefore 481 ! we have to compute the union of the k-points. 482 nkcalc = 0 483 do ii=1,tmp_nkpt 484 found = .False. 485 do spin=1,nsppol 486 if (ib_work(1, ii, spin) <= ib_work(2, ii, spin)) then 487 found = .True.; exit 488 end if 489 end do 490 if (found) then 491 nkcalc = nkcalc + 1 492 kpos(nkcalc) = ii 493 end if 494 end do 495 496 ABI_MALLOC(kcalc, (3, nkcalc)) 497 ABI_MALLOC(bstart_ks, (nkcalc, nsppol)) 498 ABI_MALLOC(nbcalc_ks, (nkcalc, nsppol)) 499 500 do ic=1,nkcalc 501 ! Index in the ib_work array 502 ii = kpos(ic) 503 ! Index in ebands. 504 ik = sigmak2ebands(ii) 505 kcalc(:,ic) = ebands%kptns(:,ik) 506 do spin=1,nsppol 507 bstart_ks(ic, spin) = 0 508 nbcalc_ks(ic, spin) = 0 509 if (ib_work(1, ii, spin) <= ib_work(2, ii, spin)) then 510 bstart_ks(ic, spin) = ib_work(1, ii, spin) 511 nbcalc_ks(ic, spin) = ib_work(2, ii, spin) - ib_work(1, ii, spin) + 1 512 end if 513 if (nbcalc_ks(ic, spin) == 0) then 514 ABI_WARNING("Spin-polarized case with nbcalc_ks == 0, don't know if code can handle it!") 515 end if 516 end do 517 end do 518 519 if (my_rank == master) then 520 ! Write info about k-points used in the calculation. 521 write(msg, "(a, i0, a, 2(f6.3, 1x), a)") & 522 " Found ", nkcalc, " k-points within sigma_erange: ", dtset%sigma_erange(:) * Ha_eV, " (eV)" 523 call wrtout(unts, msg) 524 if (any(dtset%sigma_ngkpt /= 0)) then 525 call wrtout(unts, sjoin(" These k-points belong to the sigma_ngkpt k-mesh:", ltoa(dtset%sigma_ngkpt))) 526 end if 527 write(msg, "(2(a, i0))")" min(nbcalc_ks): ", minval(nbcalc_ks), " Max(nbcalc_ks): ", maxval(nbcalc_ks) 528 call wrtout(unts, msg) 529 end if 530 531 ABI_FREE(ib_work) 532 ABI_FREE(sigmak2ebands) 533 534 end subroutine sigtk_kcalc_from_erange
m_sigtk/sigtk_kcalc_from_gaps [ Functions ]
[ Top ] [ m_sigtk ] [ Functions ]
NAME
sigtk_kcalc_from_gaps
FUNCTION
Select list of k-points and bands for self-energy matrix elements so that fundamental and direct gaps are included. INPUT dtset<dataset_type>=All input variables for this dataset. ebands<ebands_t>=The GS KS band structure (energies, occupancies, k-weights...) gaps<gaps_t>=Store location of gaps.
OUTPUT
nkcalc: Number of k-points in self-energy matrix elements. kcalc(3, nkcalc): List of k-points where the self-energy is computed. bstart_ks(nkcalc, nsppol): Initial KS band index included in self-energy matrix elements for each k-point in kcalc. nbcalc_ks(nkcalc, nsppol): Number of bands included in self-energy matrix elements for each k-point in kcalc.
SOURCE
263 subroutine sigtk_kcalc_from_gaps(dtset, ebands, gaps, nkcalc, kcalc, bstart_ks, nbcalc_ks) 264 265 !Arguments ------------------------------------ 266 type(dataset_type),intent(in) :: dtset 267 type(ebands_t),intent(in) :: ebands 268 type(gaps_t) :: gaps 269 integer,intent(out) :: nkcalc 270 !arrays 271 real(dp),allocatable,intent(out) :: kcalc(:,:) 272 integer,allocatable,intent(out) :: bstart_ks(:,:) 273 integer,allocatable,intent(out) :: nbcalc_ks(:,:) 274 275 !Local variables ------------------------------ 276 !scalars 277 integer :: spin, nsppol, ii, ik, nk_found, ifo, jj 278 logical :: found 279 !arrays 280 integer :: val_indices(ebands%nkpt, ebands%nsppol) 281 integer :: kpos(6) 282 283 ! ************************************************************************* 284 285 ABI_UNUSED((/dtset%natom/)) 286 287 call wrtout(std_out, " Including direct and fundamental KS gap in Sigma_nk") 288 ABI_CHECK(maxval(gaps%ierr) == 0, "qprange 0 cannot be used because I cannot find the gap (gap_err !=0)") 289 290 nsppol = ebands%nsppol 291 val_indices = ebands_get_valence_idx(ebands) 292 293 ! Include the direct and the fundamental KS gap. 294 ! The problem here is that kptgw and nkptgw do not depend on the spin and therefore 295 ! we have compute the union of the k-points where the fundamental and the direct gaps are located. 296 nk_found = 1; kpos(1) = gaps%fo_kpos(1,1) 297 298 ! Find the list of `interesting` kpoints. 299 do spin=1,nsppol 300 do ifo=1,3 301 ik = gaps%fo_kpos(ifo, spin) 302 found = .False.; jj = 0 303 do while (.not. found .and. jj < nk_found) 304 jj = jj + 1; found = (kpos(jj) == ik) 305 end do 306 if (.not. found) then 307 nk_found = nk_found + 1; kpos(nk_found) = ik 308 end if 309 end do 310 end do 311 312 ! Now we can define the list of k-points and the bands range. 313 nkcalc = nk_found 314 ABI_MALLOC(kcalc, (3, nkcalc)) 315 ABI_MALLOC(bstart_ks, (nkcalc, nsppol)) 316 ABI_MALLOC(nbcalc_ks, (nkcalc, nsppol)) 317 318 do ii=1,nkcalc 319 ik = kpos(ii) 320 kcalc(:,ii) = ebands%kptns(:,ik) 321 do spin=1,nsppol 322 bstart_ks(ii,spin) = val_indices(ik,spin) 323 nbcalc_ks(ii,spin) = 2 324 end do 325 end do 326 327 end subroutine sigtk_kcalc_from_gaps
m_sigtk/sigtk_kcalc_from_nkptgw [ Functions ]
[ Top ] [ m_sigtk ] [ Functions ]
NAME
sigtk_kcalc_from_nkptgw
FUNCTION
Initialize list of k-points and bands for self-energy matrix elements from nkptgw. INPUT dtset<dataset_type>=All input variables for this dataset. mband: Max number of bands.
OUTPUT
nkcalc: Number of k-points in self-energy matrix elements. kcalc(3, nkcalc): List of k-points where the self-energy is computed. bstart_ks(nkcalc, nsppol): Initial KS band index included in self-energy matrix elements for each k-point in kcalc. nbcalc_ks(nkcalc, nsppol): Number of bands included in self-energy matrix elements for each k-point in kcalc.
SOURCE
95 subroutine sigtk_kcalc_from_nkptgw(dtset, mband, nkcalc, kcalc, bstart_ks, nbcalc_ks) 96 97 !Arguments ------------------------------------ 98 type(dataset_type),intent(in) :: dtset 99 integer,intent(in) :: mband 100 integer,intent(out) :: nkcalc 101 !arrays 102 real(dp),allocatable,intent(out) :: kcalc(:,:) 103 integer,allocatable,intent(out) :: bstart_ks(:,:) 104 integer,allocatable,intent(out) :: nbcalc_ks(:,:) 105 106 !Local variables ------------------------------ 107 !scalars 108 integer :: spin, ierr, ikcalc 109 character(len=500) :: msg 110 111 ! ************************************************************************* 112 113 call wrtout(std_out, " Generating list of k-points for self-energy from kptgw and bdgw.") 114 115 nkcalc = dtset%nkptgw 116 ABI_MALLOC(kcalc, (3, nkcalc)) 117 ABI_MALLOC(bstart_ks, (nkcalc, dtset%nsppol)) 118 ABI_MALLOC(nbcalc_ks, (nkcalc, dtset%nsppol)) 119 120 kcalc = dtset%kptgw(:,1:nkcalc) 121 do spin=1,dtset%nsppol 122 bstart_ks(:,spin) = dtset%bdgw(1,1:nkcalc,spin) 123 nbcalc_ks(:,spin) = dtset%bdgw(2,1:nkcalc,spin) - dtset%bdgw(1,1:nkcalc,spin) + 1 124 end do 125 126 ! Consistency check on bdgw and mband 127 ierr = 0 128 do spin=1,dtset%nsppol 129 do ikcalc=1,nkcalc 130 if (dtset%bdgw(2,ikcalc,spin) > mband) then 131 ierr = ierr + 1 132 write(msg,'(a,2(i0,1x),2(a,i0))')& 133 "For (k, s) ",ikcalc,spin," bdgw= ",dtset%bdgw(2,ikcalc,spin), " > mband = ",mband 134 ABI_WARNING(msg) 135 end if 136 end do 137 end do 138 ABI_CHECK(ierr == 0, "Not enough bands in WFK file. See messages above. Aborting now.") 139 140 end subroutine sigtk_kcalc_from_nkptgw
m_sigtk/sigtk_kcalc_from_qprange [ Functions ]
[ Top ] [ m_sigtk ] [ Functions ]
NAME
sigtk_kcalc_from_qprange
FUNCTION
Use qprange to select the interesting k-points and the corresponding bands.
0 --> Compute the QP corrections only for the fundamental and the direct gap.
+num --> Compute the QP corrections for all the k-points in the irreducible zone and include `num`
bands above and below the Fermi level.
-num --> Compute the QP corrections for all the k-points in the irreducible zone.
Include all occupied states and `num` empty states.
INPUT
dtset<dataset_type>=All input variables for this dataset.
cryst<crystal_t>=Crystalline structure
ebands<ebands_t>=The GS KS band structure (energies, occupancies, k-weights...)
qprange: See above description
OUTPUT
nkcalc: Number of k-points in self-energy matrix elements. kcalc(3, nkcalc): List of k-points where the self-energy is computed. bstart_ks(nkcalc, nsppol): Initial KS band index included in self-energy matrix elements for each k-point in kcalc. nbcalc_ks(nkcalc, nsppol): Number of bands included in self-energy matrix elements for each k-point in kcalc.
SOURCE
170 subroutine sigtk_kcalc_from_qprange(dtset, cryst, ebands, qprange, nkcalc, kcalc, bstart_ks, nbcalc_ks) 171 172 !Arguments ------------------------------------ 173 type(dataset_type),intent(in) :: dtset 174 type(crystal_t),intent(in) :: cryst 175 type(ebands_t),intent(in) :: ebands 176 integer,intent(in) :: qprange 177 integer,intent(out) :: nkcalc 178 !arrays 179 real(dp),allocatable,intent(out) :: kcalc(:,:) 180 integer,allocatable,intent(out) :: bstart_ks(:,:) 181 integer,allocatable,intent(out) :: nbcalc_ks(:,:) 182 183 !Local variables ------------------------------ 184 !scalars 185 integer :: spin, ik, bstop, mband, sigma_nkbz 186 !arrays 187 integer :: kptrlatt(3,3) 188 integer :: val_indices(ebands%nkpt, ebands%nsppol) 189 real(dp),allocatable :: sigma_wtk(:),sigma_kbz(:,:) 190 191 ! ************************************************************************* 192 193 mband = ebands%mband 194 195 val_indices = ebands_get_valence_idx(ebands) 196 197 if (any(dtset%sigma_ngkpt /= 0)) then 198 call wrtout(std_out, " Generating list of k-points for self-energy from sigma_ngkpt and qprange.") 199 ABI_CHECK(qprange /= 0, "qprange must be != 0") 200 ! Get %kcalc from sigma_ngkpt 201 kptrlatt = 0 202 kptrlatt(1,1) = dtset%sigma_ngkpt(1); kptrlatt(2,2) = dtset%sigma_ngkpt(2); kptrlatt(3,3) = dtset%sigma_ngkpt(3) 203 call kpts_ibz_from_kptrlatt(cryst, kptrlatt, dtset%kptopt, dtset%sigma_nshiftk, dtset%sigma_shiftk, & 204 nkcalc, kcalc, sigma_wtk, sigma_nkbz, sigma_kbz) 205 ABI_FREE(sigma_kbz) 206 ABI_FREE(sigma_wtk) 207 else 208 ! Include all the k-points in the IBZ. 209 ! Note that kcalc == ebands%kptns so we can use a single ik index in the loop over k-points. 210 ! No need to map kcalc onto ebands%kptns. 211 call wrtout(std_out, " nkptgw set to 0 ==> Include all k-points in the IBZ for Sigma_nk.") 212 nkcalc = ebands%nkpt 213 ABI_MALLOC(kcalc, (3, nkcalc)) 214 kcalc = ebands%kptns 215 end if 216 217 ABI_MALLOC(bstart_ks, (nkcalc, dtset%nsppol)) 218 ABI_MALLOC(nbcalc_ks, (nkcalc, dtset%nsppol)) 219 220 if (qprange > 0) then 221 call wrtout(std_out, " Using buffer of bands above and below the Fermi level.") 222 do spin=1,dtset%nsppol 223 do ik=1,nkcalc 224 bstart_ks(ik,spin) = max(val_indices(ik,spin) - qprange, 1) 225 bstop = min(val_indices(ik,spin) + qprange, mband) 226 nbcalc_ks(ik,spin) = bstop - bstart_ks(ik,spin) + 1 227 end do 228 end do 229 230 else 231 call wrtout(std_out, " Including all occupied states and -qprange empty states.") 232 bstart_ks = 1 233 do spin=1,dtset%nsppol 234 do ik=1,nkcalc 235 nbcalc_ks(ik,spin) = min(val_indices(ik,spin) - qprange, mband) 236 end do 237 end do 238 end if 239 240 end subroutine sigtk_kcalc_from_qprange
m_sigtk/sigtk_sigma_tables [ Functions ]
[ Top ] [ m_sigtk ] [ Functions ]
NAME
sigtk_sigma_tables
FUNCTION
Build tables with the band indices used to compute the matrix elements of sigma_x and sigma_c taking into account the kind of self-energies and symmetries from esymm.
INPUTS
nkcalc: Number of k-points to compute nkibz: Number of k-points in the IBZ nsppol: Number of spins bstart_ks, bstop_ks: First and last band for each (ikcalc, spin) kcalc2ibz: Mapping kcalc --> IBZ only_diago: True if only diagonal matrix elements are wanted sigc_is_herm: True is Sigma_c is Hermitian [esymm]: Band symmetries
OUTPUT
sigxij_tab, sigcij_tab
SOURCE
840 subroutine sigtk_sigma_tables(nkcalc, nkibz, nsppol, bstart_ks, bstop_ks, kcalc2ibz, & 841 only_diago, sigc_is_herm, sigxij_tab, sigcij_tab, esymm) 842 843 use m_gwdefs, only : sigijtab_t, sigijtab_free 844 use m_esymm, only : esymm_t, esymm_failed 845 846 !Arguments ------------------------------------ 847 integer,intent(in) :: nkcalc, nkibz, nsppol 848 integer,intent(in) :: bstart_ks(nkcalc, nsppol), bstop_ks(nkcalc, nsppol) 849 logical,intent(in) :: only_diago, sigc_is_herm 850 integer,intent(in) :: kcalc2ibz(nkcalc) 851 type(sigijtab_t),allocatable,intent(inout) :: Sigxij_tab(:,:), Sigcij_tab(:,:) 852 type(esymm_t),optional,intent(in) :: esymm(nkibz, nsppol) 853 854 !Local variables------------------------------- 855 !scalars 856 integer :: spin,ikcalc,ik_ibz,bmin,bmax,bcol,brow 857 integer :: ii,idx_x,idx_c,irr_idx1,irr_idx2 858 !arrays 859 integer,allocatable :: sigc_bidx(:), sigx_bidx(:) 860 logical :: use_sym_at(nkibz, nsppol) 861 862 ! ************************************************************************* 863 864 if (allocated(Sigxij_tab)) then 865 call sigijtab_free(Sigxij_tab) 866 ABI_FREE(Sigxij_tab) 867 end if 868 if (allocated(Sigcij_tab)) then 869 call sigijtab_free(Sigcij_tab) 870 ABI_FREE(Sigcij_tab) 871 end if 872 873 ABI_MALLOC(Sigcij_tab, (nkcalc, nsppol)) 874 ABI_MALLOC(Sigxij_tab, (nkcalc, nsppol)) 875 876 use_sym_at = .FALSE. 877 if (present(esymm)) then 878 ! Create the Sig_ij tables taking advantage of the classification of the bands. 879 do spin=1,nsppol 880 do ikcalc=1,nkcalc 881 ik_ibz = kcalc2ibz(ikcalc) 882 use_sym_at(ik_ibz, spin) = .not. esymm_failed(esymm(ik_ibz, spin)) 883 end do 884 end do 885 end if 886 887 do spin=1,nsppol 888 do ikcalc=1,nkcalc 889 ik_ibz = kcalc2ibz(ikcalc) 890 891 if (use_sym_at(ik_ibz, spin)) then 892 if (only_diago) then 893 ABI_ERROR("You should not be here!") 894 end if 895 896 bmin = bstart_ks(ikcalc, spin); bmax = bstop_ks(ikcalc, spin) 897 ABI_MALLOC(Sigxij_tab(ikcalc, spin)%col, (bmin:bmax)) 898 ABI_MALLOC(Sigcij_tab(ikcalc, spin)%col, (bmin:bmax)) 899 900 do bcol=bmin,bmax 901 ABI_MALLOC(sigc_bidx, (bmax - bmin + 1)) 902 ABI_MALLOC(sigx_bidx, (bmax - bmin + 1)) 903 904 if (esymm(ik_ibz,spin)%err_status /= 0) then 905 ! Band classification failed. 906 sigc_bidx = [(ii, ii=bmin, bmax)] 907 idx_c = bmax - bmin + 1 908 sigx_bidx = [(ii,ii=bmin,bcol)] ! Hermitian 909 idx_x = bcol - bmin + 1 910 else 911 irr_idx2 = esymm(ik_ibz,spin)%b2irrep(bcol) 912 idx_c = 0 913 do brow=bmin,bmax 914 irr_idx1 = esymm(ik_ibz,spin)%b2irrep(brow) 915 if (sigc_is_herm .and. bcol < brow) CYCLE ! Only the upper triangle for HF, SEX, or COHSEX. 916 if (irr_idx1 == irr_idx2) then ! same character, add this row to the list. 917 idx_c = idx_c + 1 918 sigc_bidx(idx_c) = brow 919 end if 920 end do 921 idx_x = 0 922 do brow=bmin,bcol 923 irr_idx1 = esymm(ik_ibz,spin)%b2irrep(brow) 924 if (bcol<brow) CYCLE ! Sig_x is always Hermitian. 925 if (irr_idx1 == irr_idx2) then ! same character, add this row to the list. 926 idx_x = idx_x +1 927 sigx_bidx(idx_x) = brow 928 end if 929 end do 930 end if 931 932 ! Table for Sigma_x matrix elements taking into account symmetries of the bands. 933 ABI_MALLOC(Sigxij_tab(ikcalc, spin)%col(bcol)%bidx, (idx_x)) 934 935 Sigxij_tab(ikcalc, spin)%col(bcol)%size1 = idx_x 936 Sigxij_tab(ikcalc, spin)%col(bcol)%bidx(:) = sigx_bidx(1:idx_x) 937 !write(std_out,*)" Sigxij_tab: ikcalc, spin, bcol ",ikcalc,spin,bcol 938 !write(std_out,*)" size: ",idx_x,(Sigxij_tab(ikcalc,spin)%col(bcol)%bidx(ii),ii=1,idx_x) 939 ! 940 ! Table for Sigma_c matrix elements taking into account symmetries of the bands. 941 ABI_MALLOC(Sigcij_tab(ikcalc, spin)%col(bcol)%bidx, (idx_c)) 942 943 Sigcij_tab(ikcalc, spin)%col(bcol)%size1= idx_c 944 Sigcij_tab(ikcalc, spin)%col(bcol)%bidx(:) = sigc_bidx(1:idx_c) 945 !write(std_out,*)" Sigcij_tab: ikcalc, spin, bcol ",ikcalc,spin,bcol 946 !write(std_out,*)" size: ",idx_c,(Sigcij_tab(ikcalc,spin)%col(bcol)%bidx(ii), ii=1,idx_c) 947 948 ABI_FREE(sigx_bidx) 949 ABI_FREE(sigc_bidx) 950 end do ! bcol 951 952 else 953 ! Symmetries cannot be used for this (k,s). 954 bmin = bstart_ks(ikcalc, spin); bmax = bstop_ks(ikcalc, spin) 955 ABI_MALLOC(Sigcij_tab (ikcalc, spin)%col, (bmin:bmax)) 956 ABI_MALLOC(Sigxij_tab (ikcalc, spin)%col, (bmin:bmax)) 957 958 if (only_diago) then 959 ! QP wavefunctions == KS, therefore only diagonal elements are calculated. 960 do bcol=bmin,bmax 961 ABI_MALLOC(Sigcij_tab(ikcalc, spin)%col(bcol)%bidx, (1:1)) 962 Sigcij_tab(ikcalc, spin)%col(bcol)%size1= 1 963 Sigcij_tab(ikcalc, spin)%col(bcol)%bidx(1) = bcol 964 965 ABI_MALLOC(Sigxij_tab(ikcalc, spin)%col(bcol)%bidx, (1:1)) 966 Sigxij_tab(ikcalc, spin)%col(bcol)%size1 = 1 967 Sigxij_tab(ikcalc, spin)%col(bcol)%bidx(1) = bcol 968 end do 969 else 970 ! Use QP wavefunctions, Sigma_ij matrix is sparse but we have to classify the states in sigma. 971 ! The only thing we can do here is filling the entire matrix taking advantage of Hermiticity (if any). 972 do bcol=bmin,bmax 973 ABI_MALLOC(Sigxij_tab(ikcalc, spin)%col(bcol)%bidx, (bcol-bmin+1)) 974 Sigxij_tab(ikcalc, spin)%col(bcol)%size1= bcol-bmin+1 975 Sigxij_tab(ikcalc, spin)%col(bcol)%bidx(:) = [(ii, ii=bmin,bcol)] ! Sigma_x is Hermitian. 976 !write(std_out,*)"Sigxij_tab: ikcalc, spin, bcol ",ikcalc,spin,bcol,Sigxij_tab(ikcalc,spin)%col(bcol)%bidx(:) 977 978 ABI_MALLOC(sigc_bidx, (bmax-bmin+1)) 979 idx_c = 0 980 do brow=bmin,bmax 981 if (sigc_is_herm .and. bcol < brow) CYCLE ! Only the upper triangle of Sigc_ij is needed (SEX, COHSEX). 982 idx_c = idx_c +1 983 sigc_bidx(idx_c) = brow 984 end do 985 ABI_MALLOC(Sigcij_tab(ikcalc, spin)%col(bcol)%bidx,(idx_c)) 986 Sigcij_tab(ikcalc, spin)%col(bcol)%size1= idx_c 987 Sigcij_tab(ikcalc, spin)%col(bcol)%bidx(:) = sigc_bidx(1:idx_c) 988 ABI_FREE(sigc_bidx) 989 !write(std_out,*)"Sigcij_tab: ikcalc, spin, bcol ",ikcalc,spin,bcol,Sigcij_tab(ikcalc,spin)%col(bcol)%bidx(:) 990 end do 991 end if 992 end if 993 994 end do !ikcalc 995 end do !spin 996 997 end subroutine sigtk_sigma_tables