TABLE OF CONTENTS
- ABINIT/m_krank
- m_krank/get_ibz2bz
- m_krank/get_rank
- m_krank/krank_copy
- m_krank/krank_free
- m_krank/krank_from_kptrlatt
- m_krank/krank_get_index
- m_krank/krank_get_mapping
- m_krank/krank_new
- m_krank/krank_print
- m_krank/krank_t
- m_krank/star_from_ibz_idx
ABINIT/m_krank [ Modules ]
NAME
m_krank
FUNCTION
This module deals with rank objects for hashing k-point vector lists
COPYRIGHT
Copyright (C) 2010-2024 ABINIT group (MVer, HM, MG) This file is distributed under the terms of the GNU General Public Licence, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt . For the initials of contributors, see ~abinit/doc/developers/contributors.txt .
SOURCE
18 #if defined HAVE_CONFIG_H 19 #include "config.h" 20 #endif 21 22 #include "abi_common.h" 23 24 ! TODO: Remove file 25 !#include "libtetra.h" 26 27 module m_krank 28 29 use defs_basis 30 use m_abicore 31 use m_errors 32 33 implicit none 34 35 private
m_krank/get_ibz2bz [ Functions ]
[ Top ] [ m_krank ] [ Functions ]
NAME
get_ibz2bz
FUNCTION
Return array with the index of the IBZ wave vectors in the BZ.
INPUTS
OUTPUT
SOURCE
653 subroutine get_ibz2bz(nibz, nbz, bz2ibz, ibz2bz, ierr) 654 655 !Arguments ------------------------------------ 656 !scalars 657 integer,intent(in) :: nibz, nbz 658 integer,intent(in) :: bz2ibz(6, nbz) 659 integer,intent(out) :: ierr 660 !arrays 661 integer,allocatable,intent(out) :: ibz2bz(:) 662 663 !Local variables------------------------------- 664 !scalars 665 integer :: iq_bz, iq_ibz, isym_q, trev_q, cnt, g0_q(3) 666 logical :: isirr_q 667 668 !---------------------------------------------------------------------- 669 670 ABI_MALLOC(ibz2bz, (nibz)) 671 672 cnt = 0 673 do iq_bz=1,nbz 674 iq_ibz = bz2ibz(1, iq_bz); isym_q = bz2ibz(2, iq_bz) 675 trev_q = bz2ibz(6, iq_bz); g0_q = bz2ibz(3:5,iq_bz) 676 isirr_q = (isym_q == 1 .and. trev_q == 0 .and. all(g0_q == 0)) 677 if (isirr_q) then 678 cnt = cnt + 1 679 ibz2bz(iq_ibz) = iq_bz 680 end if 681 end do 682 683 ierr = merge(0, 1, cnt == nibz) 684 685 end subroutine get_ibz2bz
m_krank/get_rank [ Functions ]
[ Top ] [ m_krank ] [ Functions ]
NAME
get_rank
FUNCTION
This routine calculates the rank for one kpt
INPUTS
kpt = coordinates of kpoints krank = rank object for the k-grid we are using
OUTPUT
rank = rank of the kpoint
SOURCE
291 integer function get_rank(krank, kpt) result(rank) 292 293 !Arguments ------------------------------------ 294 !scalars 295 class(krank_t), intent(in) :: krank 296 !arrays 297 real(dp),intent(in) :: kpt(3) 298 299 !Local variables------------------------------- 300 !scalars 301 character(len=500) :: msg 302 !arrays 303 real(dp) :: redkpt(3) 304 305 ! ************************************************************************* 306 307 ! wrap to [0, 1[ -> replaced call to wrap2_zero2one inline, to encapsulate this module 308 if (kpt(1)>zero) then 309 redkpt(1)=mod((kpt(1)+tol12),one)-tol12 310 else 311 redkpt(1)=-mod(-(kpt(1)-one+tol12),one)+one-tol12 312 end if 313 if(abs(redkpt(1))<tol12)redkpt(1)=zero 314 315 if (kpt(2)>zero) then 316 redkpt(2)=mod((kpt(2)+tol12),one)-tol12 317 else 318 redkpt(2)=-mod(-(kpt(2)-one+tol12),one)+one-tol12 319 end if 320 if(abs(redkpt(2))<tol12)redkpt(2)=zero 321 322 if (kpt(3)>zero) then 323 redkpt(3)=mod((kpt(3)+tol12),one)-tol12 324 else 325 redkpt(3)=-mod(-(kpt(3)-one+tol12),one)+one-tol12 326 end if 327 if(abs(redkpt(3))<tol12)redkpt(3)=zero 328 329 ! rank = int(real(krank%max_linear_density)*(redkpt(3)+half+tol8 + & 330 ! real(krank%max_linear_density)*(redkpt(2)+half+tol8 + & 331 ! real(krank%max_linear_density)*(redkpt(1)+half+tol8)))) 332 rank = nint(real(krank%max_linear_density)*(redkpt(1)+half+tol8 + & 333 real(krank%max_linear_density)*(redkpt(2)+half+tol8 + & 334 real(krank%max_linear_density)*(redkpt(3)+half+tol8)))) 335 336 if (rank > krank%max_rank) then 337 write(msg,'(2(a,i0))') ' Rank should be <= max_rank: ', krank%max_rank, ' but got: ', rank 338 ABI_ERROR(msg) 339 end if 340 if (rank < krank%min_rank) then 341 !print *, "redkpt", redkpt 342 write(msg,'(2(a,i0))') ' Rank should be >= min_rank ', krank%min_rank, ' but got: ', rank 343 ABI_ERROR(msg) 344 end if 345 346 end function get_rank
m_krank/krank_copy [ Functions ]
[ Top ] [ m_krank ] [ Functions ]
NAME
krank_copy
FUNCTION
Deep copy of the object
INPUTS
OUTPUT
krank = object containing ranking and inverse ranking, to be deallocated
SOURCE
404 type(krank_t) function krank_copy(krank_in) result(krank_out) 405 406 !Arguments ------------------------------------ 407 !scalars 408 class(krank_t), intent(in) :: krank_in 409 410 ! ********************************************************************* 411 412 krank_out%max_linear_density = krank_in%max_linear_density 413 krank_out%min_rank = krank_in%min_rank 414 krank_out%max_rank = krank_in%max_rank 415 krank_out%npoints = krank_in%npoints 416 417 ABI_MALLOC(krank_out%invrank, (krank_out%min_rank:krank_out%max_rank)) 418 krank_out%invrank = krank_in%invrank 419 420 ! This is why I call it deep copy! 421 ABI_MALLOC(krank_out%kpts, (3, size(krank_in%kpts, dim=2))) 422 krank_out%kpts = krank_in%kpts 423 krank_out%kpts_owns_memory = .True. 424 425 end function krank_copy
m_krank/krank_free [ Functions ]
[ Top ] [ m_krank ] [ Functions ]
NAME
krank_free
FUNCTION
This routine deallocates the arrays in a krank_t structure
INPUTS
krank = object containing ranking and inverse ranking, to be deallocated
SOURCE
442 subroutine krank_free(krank) 443 444 !Arguments ------------------------------------ 445 class(krank_t), intent(inout) :: krank 446 447 ! ********************************************************************* 448 449 ABI_SFREE(krank%invrank) 450 ABI_SFREE(krank%rank2ikpt_) 451 ABI_SFREE(krank%rank2symtime_) 452 453 if (krank%kpts_owns_memory) then 454 ABI_SFREE_PTR(krank%kpts) 455 krank%kpts_owns_memory = .False. 456 krank%kpts => null() 457 else 458 krank%kpts => null() 459 end if 460 461 end subroutine krank_free
m_krank/krank_from_kptrlatt [ Functions ]
[ Top ] [ m_krank ] [ Functions ]
NAME
krank_from_kptrlatt
FUNCTION
This routine sets up the kpt ranks for comparing kpts
INPUTS
npt = number of kpoints (eventually irreducible) kpt = coordinates of kpoints
OUTPUT
krank = object containing ranking and inverse ranking
NOTES
By default, the object holds a reference to kpts so do not change/deallocate this array while using krank.
SOURCE
123 type(krank_t) function krank_from_kptrlatt(nkpt, kpts, kptrlatt, compute_invrank) result(new) 124 125 !Arguments ------------------------------------ 126 !scalars 127 integer,intent(in) :: nkpt 128 logical,optional,intent(in) :: compute_invrank 129 !arrays 130 integer,intent(in) :: kptrlatt(3,3) 131 real(dp),target,intent(in) :: kpts(3,nkpt) 132 133 !Local variables ------------------------- 134 !scalars 135 integer :: ii, jj, max_linear_density 136 logical :: compute_invrank_ 137 138 ! ********************************************************************* 139 140 do jj=1,3 141 do ii=1,3 142 if (ii == jj .and. kptrlatt(ii, ii) == 0) then 143 ABI_ERROR("kptrlatt with zero matrix element on the diagonal!") 144 end if 145 if (ii /= jj .and. kptrlatt(ii, jj) /= 0) then 146 ABI_ERROR("kptrlatt with non-zero off-diagonal matrix elements is not supported") 147 end if 148 end do 149 end do 150 151 compute_invrank_ = .True.; if (present(compute_invrank)) compute_invrank_ = compute_invrank 152 153 max_linear_density = maxval([kptrlatt(1,1), kptrlatt(2,2), kptrlatt(3,3)]) 154 new = krank_new(nkpt, kpts, max_linear_density=max_linear_density, compute_invrank=compute_invrank_) 155 156 end function krank_from_kptrlatt
m_krank/krank_get_index [ Functions ]
[ Top ] [ m_krank ] [ Functions ]
NAME
krank_get_index
FUNCTION
Return the index of the k-point `kpt` in the initial set. -1 if not found.
INPUTS
krank = rank object for the k-grid we are using kpt = coordinates of kpoints
OUTPUT
SOURCE
366 integer function krank_get_index(krank, kpt) result(ikpt) 367 368 !Arguments ------------------------------------ 369 !scalars 370 class(krank_t), intent(in) :: krank 371 !arrays 372 real(dp),intent(in) :: kpt(3) 373 374 !Local variables------------------------------- 375 !scalars 376 integer :: kpt_rank 377 378 ! ************************************************************************* 379 380 kpt_rank = krank%get_rank(kpt) 381 ikpt = -1 382 if (kpt_rank < krank%max_rank) ikpt = krank%invrank(kpt_rank) 383 384 end function krank_get_index
m_krank/krank_get_mapping [ Functions ]
[ Top ] [ m_krank ] [ Functions ]
NAME
krank_get_mapping
FUNCTION
Use symmetries to map input kptn2 to the list of k-points used to generate krank_t. Similar to listkk but, unlike listkk, this algo does not try to minimize the distance Mainly used to map two set of k-points associated to the same grid (e.g. BZ --> IBZ, IBZ(q) --> IBZ etc. Must faster than listkk for dense meshes although this routine requires the allocation of temporary array of shape (2, self%min_rank:self%max_rank) Returns indirect indexing list indkk.
INPUTS
kptns2(3,nkpt2)=list of final k points nkpt2=number of final k points gmet(3,3)=reciprocal space metric (bohr^-2) nsym=number of symmetry elements symafm(nsym)=(anti)ferromagnetic part of symmetry operations symmat(3,3,nsym)=symmetry operations (symrel or symrec, depending on value of use_symrec) timrev=1 if the use of time-reversal is allowed; 0 otherwise [use_symrec]: if present and true, symmat assumed to be symrec, otherwise assumed to be symrel (default) [qpt]: q-point to be added to kptns2. 0 if not specified.
OUTPUT
dksqmax=maximal value of the norm**2 of the difference between a kpt2 vector and the closest k-point found from the kptns1 set, using symmetries. indkk(6, nkpt2)= indkk(:,1)=k point index of kptns1 indkk(:,2)=symmetry operation to be applied to kpt1, to give kpt1a indkk(:,3:5)=shift in reciprocal space to be given to kpt1a, to give kpt1b, that is the closest to kpt2. indkk(:,6)=1 if time-reversal was used to generate kpt1a from kpt1, 0 otherwise
SOURCE
540 subroutine krank_get_mapping(self, nkpt2, kptns2, dksqmax, gmet, indkk, nsym, symafm, symmat, timrev, & 541 use_symrec, qpt) ! optional 542 543 !Arguments ------------------------------------ 544 !scalars 545 class(krank_t),intent(inout) :: self 546 integer,intent(in) :: nkpt2, nsym, timrev 547 real(dp),intent(out) :: dksqmax 548 logical,optional,intent(in) :: use_symrec 549 !arrays 550 integer,intent(in) :: symafm(nsym), symmat(3,3,nsym) 551 integer,intent(out) :: indkk(6, nkpt2) 552 real(dp),intent(in) :: gmet(3,3), kptns2(3,nkpt2) 553 real(dp),optional,intent(in) :: qpt(3) 554 555 !Local variables------------------------------- 556 !scalars 557 integer :: irank, ikpt1, ikpt2, itimrev, isym, ii 558 logical :: my_use_symrec 559 !arrays 560 integer :: dkint(3), my_symmat(3, 3, nsym) 561 !integer,allocatable :: rank2ikpt(:), rank2symtime(:) 562 real(dp) :: kpt1a(3), dk(3), my_qpt(3) 563 564 ! ************************************************************************* 565 566 my_use_symrec = .False.; if (present(use_symrec)) my_use_symrec = use_symrec 567 my_qpt = zero; if (present(qpt)) my_qpt = qpt 568 569 if (my_use_symrec) then 570 ! Symrec k 571 my_symmat = symmat 572 else 573 ! Symrel^T k 574 do isym=1,nsym 575 my_symmat(:,:,isym) = transpose(symmat(:,:,isym)) 576 end do 577 end if 578 579 ABI_MALLOC_IFNOT(self%rank2symtime_, (self%min_rank:self%max_rank)) 580 ABI_MALLOC_IFNOT(self%rank2ikpt_, (self%min_rank:self%max_rank)) 581 self%rank2ikpt_ = -1 !; self%rank2symtime_ = -1 582 583 do ikpt1=1,self%npoints 584 585 do itimrev=0,timrev 586 do isym=1,nsym 587 ! Do not use magnetic symmetries. 588 if (symafm(isym) == -1) cycle 589 590 kpt1a = (1 - 2*itimrev) * matmul(my_symmat(:, :, isym), self%kpts(:, ikpt1)) 591 irank = self%get_rank(kpt1a) 592 if (self%rank2ikpt_(irank) == -1) then 593 self%rank2ikpt_(irank) = ikpt1 594 self%rank2symtime_(irank) = isym + itimrev * nsym 595 end if 596 end do 597 end do 598 599 end do 600 601 dksqmax = zero 602 do ikpt2=1,nkpt2 603 irank = self%get_rank(kptns2(:, ikpt2) + my_qpt) 604 ikpt1 = self%rank2ikpt_(irank) 605 ii = self%rank2symtime_(irank) 606 isym = 1 + mod(ii - 1, nsym) 607 itimrev = (ii - 1) / nsym 608 kpt1a = (1 - 2 * itimrev) * matmul(my_symmat(:, :, isym), self%kpts(:, ikpt1)) 609 dk(:) = kptns2(:,ikpt2) + my_qpt - kpt1a(:) 610 dkint(:) = nint(dk(:) + tol12) 611 612 indkk(1, ikpt2) = ikpt1 613 indkk(2, ikpt2) = isym 614 indkk(3:5, ikpt2) = dkint(:) 615 indkk(6, ikpt2) = itimrev 616 617 ! Compute norm of the difference vector. 618 dk(:) = dk(:) - dkint(:) 619 620 dksqmax = max(dksqmax, & 621 gmet(1,1)*dk(1)**2 + gmet(2,2)*dk(2)**2 + gmet(3,3)*dk(3)**2 + & 622 two * (gmet(2,1)*dk(2)*dk(1) + gmet(3,2)*dk(3)*dk(2)+gmet(3,1)*dk(3)*dk(1))) 623 624 !if (dksqmax > tol8) then 625 ! print *, "kbase:", self%kpts(:, ikpt1) 626 ! print *, "k", kptns2(:, ikpt2) 627 ! print *, "k + q", kptns2(:, ikpt2) + my_qpt 628 ! print *, "dk", dk, "dksqmax", dksqmax 629 !end if 630 end do 631 632 !ABI_FREE(self%rank2ikpt_) 633 !ABI_FREE(self%rank2symtime_) 634 635 end subroutine krank_get_mapping
m_krank/krank_new [ Functions ]
[ Top ] [ m_krank ] [ Functions ]
NAME
krank_new
FUNCTION
This routine sets up the kpt ranks for comparing kpts
INPUTS
npt = number of kpoints (eventually irreducible) kpt = coordinates of kpoints time_reversal = true or false to use time reversal symmetry. Default is true, but only important if nsym and symrec are present
SOURCE
174 type(krank_t) function krank_new(nkpt, kpts, nsym, symrec, time_reversal, max_linear_density, compute_invrank) result(new) 175 176 !Arguments ------------------------------------ 177 !scalars 178 integer,intent(in) :: nkpt 179 integer,intent(in), optional :: nsym 180 logical,intent(in), optional :: time_reversal 181 integer,optional,intent(in) :: max_linear_density 182 logical,optional,intent(in) :: compute_invrank 183 !arrays 184 real(dp),target,intent(in) :: kpts(3,nkpt) 185 integer,intent(in), optional :: symrec(3,3, *) 186 187 !Local variables ------------------------- 188 !scalars 189 integer :: ikpt, isym, symkptrank, irank, timrev, itim 190 logical :: compute_invrank_ 191 real(dp) :: smallestlen 192 character(len=500) :: msg 193 !arrays 194 real(dp) :: symkpt(3) 195 196 ! ********************************************************************* 197 198 compute_invrank_ = .True.; if (present(compute_invrank)) compute_invrank_ = compute_invrank 199 200 new%kpts => kpts 201 new%kpts_owns_memory = .False. 202 203 if (.not. present(max_linear_density)) then 204 ! Find smallest linear length from input kpts 205 smallestlen = one 206 do ikpt=1, nkpt 207 if (abs(kpts(1,ikpt)) > tol10) smallestlen = min(smallestlen, abs(kpts(1,ikpt))) 208 if (abs(kpts(2,ikpt)) > tol10) smallestlen = min(smallestlen, abs(kpts(2,ikpt))) 209 if (abs(kpts(3,ikpt)) > tol10) smallestlen = min(smallestlen, abs(kpts(3,ikpt))) 210 end do 211 new%max_linear_density = nint(one/smallestlen) 212 else 213 ! Get it from input 214 new%max_linear_density = max_linear_density 215 end if 216 217 new%npoints = nkpt 218 new%min_rank = nint(real(new%max_linear_density)*(half+tol8 +& 219 real(new%max_linear_density)*(half+tol8 +& 220 real(new%max_linear_density)*(half+tol8)))) 221 222 new%max_rank = nint(real(new%max_linear_density)*(1+half+tol8 +& 223 real(new%max_linear_density)*(1+half+tol8 +& 224 real(new%max_linear_density)*(1+half+tol8)))) 225 226 timrev = 2 227 new%time_reversal = .true. 228 if (present(time_reversal)) then 229 if (.not. time_reversal) timrev = 1 230 new%time_reversal = .false. 231 end if 232 233 ! Ensure kpt(i)+one is positive, and the smallest difference between kpts should be larger than 1/100 ie ngkpt < 100. 234 ! the following fills invrank for the k-points in the list provided (may be only the irred kpts) 235 if (compute_invrank_) then 236 ABI_MALLOC(new%invrank, (new%min_rank:new%max_rank)) 237 new%invrank(:) = -1 238 239 do ikpt=1,nkpt 240 irank = new%get_rank(kpts(:,ikpt)) 241 242 if (irank > new%max_rank .or. irank < new%min_rank) then 243 write(msg,'(a,2i0)')" rank above max_rank or below min_rank, ikpt, rank ", ikpt, irank 244 ABI_ERROR(msg) 245 end if 246 new%invrank(irank) = ikpt 247 end do 248 end if 249 250 ! if symrec is provided, fill invrank with appropriate irred kpt indices 251 ! for symmetry completion: kptrank_t%invrank points to the irred k-point 252 ! equivalent to the k-point whose rank is provided 253 if (present(symrec)) then 254 if(.not. present(nsym)) then 255 ABI_ERROR("need both symrec and nsym arguments together") 256 end if 257 do ikpt=1,nkpt 258 ! itim == 1 for positive, and itim==2 gives Kramers opposite of k-point 259 ! favor the former by looping it last 260 do itim = timrev, 1, -1 261 do isym = 1, nsym 262 symkpt = (-1)**(itim+1) * matmul(symrec(:,:,isym), kpts(:, ikpt)) 263 symkptrank = new%get_rank(symkpt(:)) 264 new%invrank(symkptrank) = ikpt 265 end do 266 end do 267 end do 268 end if 269 270 end function krank_new
m_krank/krank_print [ Functions ]
[ Top ] [ m_krank ] [ Functions ]
NAME
krank_print
FUNCTION
This routine prints the arrays and dimensions of a krank_t structure
INPUTS
krank = object containing ranking and inverse ranking unout = unit for open file to print to
SOURCE
480 subroutine krank_print(krank, unout) 481 482 !Arguments ------------------------------------ 483 !scalars 484 integer, intent(in) :: unout 485 !arrays 486 class(krank_t), intent(in) :: krank 487 488 ! ********************************************************************* 489 490 write(unout, *) 491 write(unout, '(a)') ' Dump of the contents of a krank_t structure with k-point rank information' 492 write(unout, '(a,i0)') ' max linear density of points in 3 directions: max_linear_density = ', krank%max_linear_density 493 write(unout, '(a,i0)') ' maximum rank for any point in grid: max_rank = ', krank%max_rank 494 write(unout, '(a,i0)') ' number of points in input grid: npoints = ', krank%npoints 495 write(unout, *) 496 write(unout, '(a)') ' invrank array = ' 497 write(unout, '(i0)') krank%invrank(:) 498 write(unout, *) 499 500 end subroutine krank_print
m_krank/krank_t [ Types ]
NAME
krank_t
FUNCTION
structure to contain a rank/inverse rank pair of arrays, with dimensions
SOURCE
47 type,public :: krank_t 48 49 integer :: max_linear_density 50 51 integer :: min_rank 52 53 integer :: max_rank 54 55 integer :: npoints 56 57 logical :: time_reversal 58 59 logical :: kpts_owns_memory = .False. 60 61 integer,allocatable :: invrank(:) 62 63 real(dp),ABI_CONTIGUOUS pointer :: kpts(:,:) 64 ! Reference to input k-points or copy of the array depending on kpts_owns_memory 65 66 ! Internal tables used by krank_get_mapping 67 integer,allocatable :: rank2ikpt_(:), rank2symtime_(:) 68 69 contains 70 71 procedure :: get_rank 72 ! Calculates the rank for one kpt 73 74 procedure :: get_index => krank_get_index 75 ! Return the index of the k-point `kpt` in the initial set. -1 if not found. 76 77 procedure :: copy => krank_copy 78 ! Deep copy of the object 79 80 procedure :: free => krank_free 81 ! Free memory 82 83 procedure :: print => krank_print 84 ! Prints the arrays and dimensions of a krank_t structure 85 86 procedure :: get_mapping => krank_get_mapping 87 ! Use symmetries to map input kptn2 to the list of k-points used to generate krank_t. 88 ! Similar to listkk but, unlike listkk, this algo does not try to minimize the distance. 89 ! Must faster than listkk for dense meshes 90 91 end type krank_t 92 93 public :: krank_from_kptrlatt ! Initialize object from kptrlatt 94 public :: krank_new ! Sets up the kpt ranks for comparing kpts 95 96 public :: get_ibz2bz ! Return array with the index of the IBZ wave vectors in the BZ. 97 public :: star_from_ibz_idx ! Return array with the indices of the star of the ik_ibz wavevector in the IBZ.
m_krank/star_from_ibz_idx [ Functions ]
[ Top ] [ m_krank ] [ Functions ]
NAME
star_from_ibz_idx
FUNCTION
Return array with the indices of the star of the ik_ibz wavevector in the IBZ. Return number of points in the star, allocate array with the indices of the star points in the BZ.
INPUTS
OUTPUT
SOURCE
705 subroutine star_from_ibz_idx(ik_ibz, nkbz, bz2ibz, nk_in_star, kstar_bz_inds) 706 707 !Arguments ------------------------------------ 708 integer,intent(in) :: ik_ibz, nkbz 709 integer,intent(out) :: nk_in_star 710 integer,intent(in) :: bz2ibz(6, nkbz) 711 integer,allocatable,intent(out) :: kstar_bz_inds(:) 712 713 !Local variables------------------------------- 714 !scalars 715 integer :: iq_bz 716 717 !---------------------------------------------------------------------- 718 719 nk_in_star = count(bz2ibz(1, :) == ik_ibz) 720 ABI_MALLOC(kstar_bz_inds, (nk_in_star)) 721 722 nk_in_star = 0 723 do iq_bz=1,nkbz 724 if (bz2ibz(1, iq_bz) /= ik_ibz) cycle 725 nk_in_star = nk_in_star + 1 726 kstar_bz_inds(nk_in_star) = iq_bz 727 end do 728 729 end subroutine star_from_ibz_idx