TABLE OF CONTENTS
ABINIT/get_full_kgrid [ Functions ]
NAME
get_full_kgrid
FUNCTION
create full grid of kpoints and find equivalent irred ones. Duplicates work in getkgrid, but need all outputs of kpt_fullbz, and indkpt
COPYRIGHT
Copyright (C) 2002-2018 ABINIT group (MVer,XG) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt . For the initials of contributors, see ~abinit/doc/developers/contributors.txt .
INPUTS
kpt(3,nkpt)=irreducible kpoints kptrlatt(3,3)=lattice vectors for full kpoint grid nkpt=number of irreducible kpoints nkpt_fullbz=number of kpoints in full brillouin zone nshiftk=number of kpoint grid shifts nsym=number of symmetries shiftk(3,nshiftk)=kpoint shifts symrel(3,3,nsym)=symmetry matrices in real space
OUTPUT
indkpt(nkpt_fullbz)=non-symmetrized indices of the k-points (see symkpt.f) kpt_fullbz(3,nkpt_fullbz)=kpoints in full brillouin zone
NOTES
MG: The present inplementation always assumes kptopt==1 !!!! TODO: This routine should be removed
PARENTS
m_phonons
CHILDREN
destroy_kptrank,get_kpt_fullbz,get_rank_1kpt,mati3inv,mkkptrank
SOURCE
45 #if defined HAVE_CONFIG_H 46 #include "config.h" 47 #endif 48 49 #include "abi_common.h" 50 51 subroutine get_full_kgrid(indkpt,kpt,kpt_fullbz,kptrlatt,nkpt,& 52 & nkpt_fullbz,nshiftk,nsym,shiftk,symrel) 53 54 use defs_basis 55 use m_kptrank 56 use m_profiling_abi 57 use m_errors 58 59 use m_numeric_tools, only : wrap2_pmhalf 60 61 !This section has been created automatically by the script Abilint (TD). 62 !Do not modify the following lines by hand. 63 #undef ABI_FUNC 64 #define ABI_FUNC 'get_full_kgrid' 65 use interfaces_32_util 66 use interfaces_56_recipspace, except_this_one => get_full_kgrid 67 !End of the abilint section 68 69 implicit none 70 71 !Arguments ------------------------------------ 72 !scalars 73 integer,intent(in) :: nkpt,nkpt_fullbz,nshiftk,nsym 74 !arrays 75 integer,intent(in) :: kptrlatt(3,3),symrel(3,3,nsym) 76 integer,intent(out) :: indkpt(nkpt_fullbz) 77 real(dp),intent(in) :: kpt(3,nkpt),shiftk(3,nshiftk) 78 real(dp),intent(out) :: kpt_fullbz(3,nkpt_fullbz) 79 80 !Local variables------------------------------- 81 !scalars 82 integer :: ikpt,isym,itim,timrev 83 integer :: symrankkpt 84 character(len=500) :: message 85 type(kptrank_type) :: kptrank_t 86 87 !arrays 88 integer :: inv_symrel(3,3,nsym) 89 real(dp) :: k2(3) 90 91 ! ********************************************************************* 92 93 !Invert symrels => gives symrels for kpoints 94 95 do isym=1,nsym 96 call mati3inv (symrel(:,:,isym),inv_symrel(:,:,isym)) 97 end do 98 99 call get_kpt_fullbz(kpt_fullbz,kptrlatt,nkpt_fullbz,nshiftk,shiftk) 100 101 !make full k-point rank arrays 102 call mkkptrank (kpt,nkpt,kptrank_t) 103 104 ! 105 !find equivalence to irred kpoints in kpt 106 ! 107 indkpt(:) = 0 108 timrev=1 ! includes the time inversion symmetry 109 do ikpt=1,nkpt_fullbz 110 do isym=1,nsym 111 do itim=1,(1-2*timrev),-2 112 113 k2(:) = itim*(inv_symrel(:,1,isym)*kpt_fullbz(1,ikpt) + & 114 & inv_symrel(:,2,isym)*kpt_fullbz(2,ikpt) + & 115 & inv_symrel(:,3,isym)*kpt_fullbz(3,ikpt)) 116 117 call get_rank_1kpt (k2,symrankkpt,kptrank_t) 118 if (kptrank_t%invrank(symrankkpt) /= -1) indkpt(ikpt) = kptrank_t%invrank(symrankkpt) 119 120 end do ! loop time reversal symmetry 121 end do ! loop sym ops 122 123 if (indkpt(ikpt) == 0) then 124 write (message,'(a,i0)')' indkpt(ikpt) is still 0: no irred kpoint is equiv to ikpt ',ikpt 125 MSG_BUG(message) 126 end if 127 end do ! loop full kpts 128 129 call destroy_kptrank (kptrank_t) 130 131 end subroutine get_full_kgrid