TABLE OF CONTENTS
ABINIT/ep_setupqpt [ Functions ]
NAME
ep_setupqpt
FUNCTION
set up qpoint grid for elphon. 2 modes, either uniform grid from anaddb input nqpt or take qpt from anaddb input (explicitly listed)
COPYRIGHT
Copyright (C) 2009-2018 ABINIT group (MVer) 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 .
INPUTS
crystal>crystal_t>=data type gathering info on the crystalline structure. anaddb_dtset=dataset with input variables %qgrid_type gives type of q grid 1=uniform 2=take from input %ep_nqpt number of auxiliary qpoints %ep_qptlist list of qpoints,
OUTPUT
PARENTS
elphon
CHILDREN
getkgrid,smpbz,symkpt,wrap2_pmhalf,wrtout
NOTES
SOURCE
38 #if defined HAVE_CONFIG_H 39 #include "config.h" 40 #endif 41 42 #include "abi_common.h" 43 44 subroutine ep_setupqpt (elph_ds,crystal,anaddb_dtset,qptrlatt,timrev) 45 46 use defs_basis 47 use defs_elphon 48 use m_errors 49 use m_profiling_abi 50 51 use m_numeric_tools, only : wrap2_pmhalf 52 use m_anaddb_dataset, only : anaddb_dataset_type 53 use m_crystal, only : crystal_t 54 55 !This section has been created automatically by the script Abilint (TD). 56 !Do not modify the following lines by hand. 57 #undef ABI_FUNC 58 #define ABI_FUNC 'ep_setupqpt' 59 use interfaces_14_hidewrite 60 use interfaces_41_geometry 61 use interfaces_56_recipspace 62 !End of the abilint section 63 64 implicit none 65 66 !Arguments ------------------------------- 67 !scalars 68 integer, intent(in) :: timrev 69 type(crystal_t),intent(in) :: crystal 70 type(anaddb_dataset_type), intent(in) :: anaddb_dtset 71 type(elph_type), intent(inout) :: elph_ds 72 !arrays 73 integer, intent(out) :: qptrlatt(3,3) 74 75 !Local variables ------------------------- 76 !scalars 77 integer :: nqshft,option,iqpt, nqpt1 78 integer :: iscf,mqpt,iout,berryopt,nqpt_computed 79 real(dp) :: qptrlen, res 80 character(len=500) :: message 81 !arrays 82 integer :: vacuum(3) 83 integer,allocatable :: indqpt1(:) 84 real(dp) :: kpt(3) 85 real(dp),allocatable :: wtq_folded(:) 86 real(dp), allocatable :: wtq(:),qpt_full(:,:),tmpshifts(:,:) 87 88 ! ********************************************************************* 89 90 !default is to expect a uniform grid 91 elph_ds%tuniformgrid = 1 92 93 !if we use the normal grid way of generating the qpoints: 94 if (anaddb_dtset%qgrid_type==1) then 95 ! qpoint lattice vectors (inverse, like kptrlatt) 96 qptrlatt(:,:)=0 97 qptrlatt(1,1)=anaddb_dtset%ngqpt(1) 98 qptrlatt(2,2)=anaddb_dtset%ngqpt(2) 99 qptrlatt(3,3)=anaddb_dtset%ngqpt(3) 100 101 if (anaddb_dtset%nqshft /= 1) then 102 ! try to reduce the qpoint grid to a single qshift, otherwise stop 103 ! dummy args for call to getkgrid 104 vacuum(:) = 0 105 iscf = 3 106 107 mqpt = anaddb_dtset%ngqpt(1)*anaddb_dtset%ngqpt(2)*anaddb_dtset%ngqpt(3)*anaddb_dtset%nqshft 108 ABI_ALLOCATE(qpt_full,(3,mqpt)) 109 ABI_ALLOCATE(wtq,(mqpt)) 110 ABI_ALLOCATE(tmpshifts,(3,210)) 111 112 wtq(:) = one 113 114 tmpshifts(:,:) = zero 115 tmpshifts(:,1:4) = anaddb_dtset%q1shft(:,:) 116 117 iout=6 118 119 berryopt = 1 120 121 ! just call with identity, to get full set of kpts in qpt_full, but 122 ! reduce qshfts 123 124 nqshft=anaddb_dtset%nqshft 125 call getkgrid(0,0,iscf,qpt_full,3,qptrlatt,qptrlen, & 126 & 1,mqpt,nqpt_computed,nqshft,1,crystal%rprimd,tmpshifts,crystal%symafm, & 127 & crystal%symrel,vacuum,wtq) 128 ABI_DEALLOCATE(qpt_full) 129 ABI_DEALLOCATE(wtq) 130 ABI_DEALLOCATE(tmpshifts) 131 132 if (anaddb_dtset%nqshft /= 1) then 133 write (message,'(a,i0)')& 134 & ' multiple qpt shifts not treated yet (should be possible), nqshft= ', anaddb_dtset%nqshft 135 MSG_ERROR(message) 136 end if 137 end if ! end multiple shifted qgrid 138 139 140 write(message,'(a,9(i0,1x))')' elphon : enter smpbz with qptrlatt = ',qptrlatt 141 call wrtout(std_out,message,'COLL') 142 143 option=1 144 ! mqpt=anaddb_dtset%ngqpt(1)*anaddb_dtset%ngqpt(2)*anaddb_dtset%ngqpt(3)*anaddb_dtset%nqshft 145 mqpt= qptrlatt(1,1)*qptrlatt(2,2)*qptrlatt(3,3) & 146 & +qptrlatt(1,2)*qptrlatt(2,3)*qptrlatt(3,1) & 147 & +qptrlatt(1,3)*qptrlatt(2,1)*qptrlatt(3,2) & 148 & -qptrlatt(1,2)*qptrlatt(2,1)*qptrlatt(3,3) & 149 & -qptrlatt(1,3)*qptrlatt(2,2)*qptrlatt(3,1) & 150 & -qptrlatt(1,1)*qptrlatt(2,3)*qptrlatt(3,2) 151 152 ABI_ALLOCATE(qpt_full,(3,mqpt)) 153 iout = 6 154 call smpbz(anaddb_dtset%brav,iout,qptrlatt,mqpt,elph_ds%nqpt_full,anaddb_dtset%nqshft,option,anaddb_dtset%q1shft,qpt_full) 155 156 157 ! save the q-grid for future reference 158 ABI_ALLOCATE(elph_ds%qpt_full,(3,elph_ds%nqpt_full)) 159 160 ! reduce qpt_full to correct zone 161 do iqpt=1,elph_ds%nqpt_full 162 call wrap2_pmhalf(qpt_full(1,iqpt),kpt(1),res) 163 call wrap2_pmhalf(qpt_full(2,iqpt),kpt(2),res) 164 call wrap2_pmhalf(qpt_full(3,iqpt),kpt(3),res) 165 qpt_full(:,iqpt) = kpt 166 elph_ds%qpt_full(:,iqpt)=kpt 167 end do 168 ABI_DEALLOCATE(qpt_full) 169 170 else if (anaddb_dtset%qgrid_type==2) then ! use explicit list of qpoints from anaddb input 171 qptrlatt(:,:)=0 172 qptrlatt(1,1)=1 173 qptrlatt(2,2)=1 174 qptrlatt(3,3)=1 175 176 elph_ds%nqpt_full=anaddb_dtset%ep_nqpt 177 ABI_ALLOCATE(elph_ds%qpt_full,(3,elph_ds%nqpt_full)) 178 179 elph_ds%qpt_full = anaddb_dtset%ep_qptlist 180 181 elph_ds%tuniformgrid = 0 182 end if ! type of qgrid for elphon 183 184 !================================================================= 185 !Calculate weights, needed to estimate lambda using the weighted 186 !sum of the uninterpolated e-ph matrix elements 187 !================================================================= 188 call wrtout(std_out,' setqgrid : calling symkpt to find irred q points',"COLL") 189 190 ABI_ALLOCATE(indqpt1,(elph_ds%nqpt_full)) 191 ABI_ALLOCATE(wtq_folded,(elph_ds%nqpt_full)) 192 ABI_ALLOCATE(wtq,(elph_ds%nqpt_full)) 193 194 wtq(:) = one/dble(elph_ds%nqpt_full) !weights normalized to unity 195 196 ! 197 !NOTE: this reduction of irred qpt may not be identical to that in GKK file 198 !which would be more practical to use. 199 ! 200 iout=0 !do not write to ab_out 201 !should we save indqpt1 for use inside elph_ds? 202 call symkpt(0,crystal%gmet,indqpt1,iout,elph_ds%qpt_full,elph_ds%nqpt_full,nqpt1,crystal%nsym,crystal%symrec,& 203 & timrev,wtq,wtq_folded) 204 205 write (message,'(2a,i0)')ch10,' Number of irreducible q-points = ',nqpt1 206 call wrtout(std_out,message,'COLL') 207 elph_ds%nqptirred=nqpt1 208 209 call wrtout(std_out,' === Irreducible q points with weights ==== ','COLL') 210 211 do iqpt=1,elph_ds%nqpt_full 212 if (wtq_folded(iqpt) /= zero) then 213 write (message,'(1x,i4,a2,4es16.8)')iqpt,') ',elph_ds%qpt_full(:,iqpt),wtq_folded(iqpt) 214 call wrtout(std_out,message,'COLL') 215 end if 216 end do 217 218 call wrtout(std_out,ch10,'COLL') 219 220 ABI_ALLOCATE(elph_ds%wtq,(elph_ds%nqpt_full)) 221 222 elph_ds%wtq(:)=wtq_folded(:) 223 !MEMO indqpt could be useful to test the qgrid read by abinit 224 ABI_DEALLOCATE(indqpt1) 225 ABI_DEALLOCATE(wtq_folded) 226 ABI_DEALLOCATE(wtq) 227 228 end subroutine ep_setupqpt