TABLE OF CONTENTS
ABINIT/integrate_gamma_tr_lova [ Functions ]
NAME
integrate_gamma_tr_lova
FUNCTION
This routine integrates the TRANSPORT electron phonon coupling matrices over the kpoints on the fermi surface. A dependency on qpoint remains for gamma_qpt_in/out Copied from integrate_gamma
COPYRIGHT
Copyright (C) 2004-2018 ABINIT group (JPC,MJV) 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
elph_ds = elphon datastructure with data and dimensions elph_ds%qpt_full = qpoint coordinates FSfullpqtofull = mapping of k+q to k
OUTPUT
elph_tr_ds%gamma_qpt_trout elph_tr_ds%gamma_qpt_trin
PARENTS
elphon
CHILDREN
wrtout,xmpi_sum
SOURCE
37 #if defined HAVE_CONFIG_H 38 #include "config.h" 39 #endif 40 41 #include "abi_common.h" 42 43 subroutine integrate_gamma_tr_lova(elph_ds,FSfullpqtofull,elph_tr_ds) 44 45 use defs_basis 46 use defs_elphon 47 use m_profiling_abi 48 use m_errors 49 use m_xmpi 50 51 !This section has been created automatically by the script Abilint (TD). 52 !Do not modify the following lines by hand. 53 #undef ABI_FUNC 54 #define ABI_FUNC 'integrate_gamma_tr_lova' 55 use interfaces_14_hidewrite 56 !End of the abilint section 57 58 implicit none 59 60 !Arguments ------------------------------------ 61 !scalars 62 type(elph_tr_type), intent(inout) :: elph_tr_ds 63 type(elph_type),intent(in) :: elph_ds 64 !arrays 65 integer,intent(in) :: FSfullpqtofull(elph_ds%k_phon%nkpt,elph_ds%nqpt_full) 66 67 !Local variables------------------------------- 68 !scalars 69 integer :: ikpt_phon,ikpt_phonq,ib1,ib2,ibeff,ierr,iqpt,iqpt_fullbz,isppol 70 integer :: itensor, icomp, jcomp,comm 71 integer :: fib1, fib2 72 integer :: ik_this_proc 73 real(dp) :: etain, etaout 74 character(len=500) :: message 75 !arrays 76 real(dp) :: elvelock(3), elvelockpq(3) 77 real(dp),allocatable :: tmp_gkk(:,:,:,:) 78 79 ! ************************************************************************* 80 81 comm = xmpi_world 82 83 ib1=elph_ds%nbranch*elph_ds%nbranch ; ib2=elph_ds%nqpt_full 84 ABI_STAT_ALLOCATE(elph_tr_ds%gamma_qpt_trin,(2,9,ib1,elph_ds%nsppol,ib2), ierr) 85 ABI_CHECK(ierr==0, 'trying to allocate array elph_tr_ds%gamma_qpt_trin') 86 elph_tr_ds%gamma_qpt_trin = zero 87 88 ABI_STAT_ALLOCATE(elph_tr_ds%gamma_qpt_trout,(2,9,ib1,elph_ds%nsppol,ib2), ierr) 89 ABI_CHECK(ierr==0, 'trying to allocate array elph_tr_ds%gamma_qpt_trout') 90 elph_tr_ds%gamma_qpt_trout = zero 91 92 !information 93 if (elph_ds%gkqwrite == 0) then 94 write (message,'(a)')' integrate_gamma_tr : keeping gamma matrices in memory' 95 call wrtout(std_out,message,'COLL') 96 else if (elph_ds%gkqwrite == 1) then 97 write (message,'(a)')' integrate_gamma_tr : reading gamma matrices from disk' 98 call wrtout(std_out,message,'COLL') 99 else 100 write (message,'(3a,i3)')' integrate_gamma_tr : BUG-',ch10,& 101 & ' Wrong value for gkqwrite = ',elph_ds%gkqwrite 102 MSG_ERROR(message) 103 end if 104 105 !allocate temp variables 106 ABI_STAT_ALLOCATE(tmp_gkk,(2,elph_ds%ngkkband**2,elph_ds%nbranch**2,elph_ds%nsppol), ierr) 107 ABI_CHECK(ierr==0, 'trying to allocate array tmp_gkkout') 108 109 do iqpt=1,elph_ds%nqptirred 110 iqpt_fullbz = elph_ds%qirredtofull(iqpt) 111 write(std_out,*)'iqpt, iqptfullbz ',iqpt, iqpt_fullbz 112 113 do ik_this_proc =1,elph_ds%k_phon%my_nkpt 114 ikpt_phon = elph_ds%k_phon%my_ikpt(ik_this_proc) 115 116 if (elph_ds%gkqwrite == 0) then 117 tmp_gkk = elph_ds%gkk_qpt(:,:,:,ik_this_proc,:,iqpt) 118 else if (elph_ds%gkqwrite == 1) then 119 read(elph_ds%unitgkq,REC=((iqpt-1)*elph_ds%k_phon%my_nkpt+ik_this_proc)) tmp_gkk 120 end if 121 122 ikpt_phonq = FSfullpqtofull(ikpt_phon,iqpt_fullbz) 123 124 do isppol=1,elph_ds%nsppol 125 do ib1=1,elph_ds%ngkkband 126 fib1=ib1+elph_ds%minFSband-1 127 elvelock(:)=elph_tr_ds%el_veloc(ikpt_phon,fib1,:,isppol) 128 129 do ib2=1,elph_ds%ngkkband 130 ibeff=ib2+(ib1-1)*elph_ds%ngkkband 131 fib2=ib2+elph_ds%minFSband-1 132 elvelockpq(:)= elph_tr_ds%el_veloc(ikpt_phonq,fib2,:,isppol) 133 134 135 ! MJV 31/03/2009: Note that the following is valid for any geometry, not just cubic! 136 ! see eq 5 and 6 of prb 36 4103 (Al-Lehaibi et al 1987) 137 ! see also Allen PRB 17 3725 138 ! generalization to tensorial quantities is simple, by keeping the directional 139 ! references of velock and velockpq as indices. 140 do icomp = 1, 3 141 do jcomp = 1, 3 142 itensor = (icomp-1)*3+jcomp 143 ! FIXME: could use symmetry i <-> j 144 145 etain = elvelock(icomp)*elvelockpq(jcomp) 146 etaout = elvelock(icomp)*elvelock(jcomp) 147 148 149 elph_tr_ds%gamma_qpt_trin(:,itensor,:,isppol,iqpt_fullbz) = & 150 & elph_tr_ds%gamma_qpt_trin(:,itensor,:,isppol,iqpt_fullbz) + & 151 & tmp_gkk(:,ibeff,:,isppol) & 152 & *etain & 153 & *elph_ds%gkk_intweight(ib1,ikpt_phon,isppol)*elph_ds%gkk_intweight(ib2,ikpt_phonq,isppol) 154 155 elph_tr_ds%gamma_qpt_trout(:,itensor,:,isppol,iqpt_fullbz) = & 156 & elph_tr_ds%gamma_qpt_trout(:,itensor,:,isppol,iqpt_fullbz) + & 157 & tmp_gkk(:,ibeff,:,isppol) & 158 & *etaout & 159 & *elph_ds%gkk_intweight(ib1,ikpt_phon,isppol)*elph_ds%gkk_intweight(ib2,ikpt_phonq,isppol) 160 161 end do 162 end do 163 end do 164 end do 165 166 end do ! isppol 167 end do ! ik 168 169 end do ! iq 170 171 ABI_DEALLOCATE(tmp_gkk) 172 173 call xmpi_sum (elph_tr_ds%gamma_qpt_trout, comm, ierr) 174 call xmpi_sum (elph_tr_ds%gamma_qpt_trin, comm, ierr) 175 176 177 ! 178 !normalize tensor with 1/sqrt(v_x**2 * v_y**2) 179 ! 180 !move the veloc into mka2f_tr_lova, where T dependence is dealt with 181 !This will cause some slight difference to the results 182 if (.true.) then 183 do isppol=1, elph_ds%nsppol 184 do icomp = 1, 3 185 do jcomp = 1, 3 186 itensor = (icomp-1)*3+jcomp 187 if(abs(elph_tr_ds%FSelecveloc_sq(icomp,isppol))>tol14**2 .and. abs(elph_tr_ds%FSelecveloc_sq(jcomp,isppol))>tol14**2)then 188 elph_tr_ds%gamma_qpt_trin(:,itensor,:,isppol,:) = elph_tr_ds%gamma_qpt_trin(:,itensor,:,isppol,:) / & 189 & sqrt(elph_tr_ds%FSelecveloc_sq(icomp,isppol)*elph_tr_ds%FSelecveloc_sq(jcomp,isppol)) 190 elph_tr_ds%gamma_qpt_trout(:,itensor,:,isppol,:) = elph_tr_ds%gamma_qpt_trout(:,itensor,:,isppol,:) / & 191 & sqrt(elph_tr_ds%FSelecveloc_sq(icomp,isppol)*elph_tr_ds%FSelecveloc_sq(jcomp,isppol)) 192 else 193 ! XG120528 Fixed problem with zero velocity 194 elph_tr_ds%gamma_qpt_trin(:,itensor,:,isppol,:)=zero 195 elph_tr_ds%gamma_qpt_trout(:,itensor,:,isppol,:)=zero 196 end if 197 end do 198 end do 199 end do ! isppol 200 end if 201 202 !need prefactor of 1/nkpt for each integration over 1 kpoint index. 203 !NOT INCLUDED IN elph_ds%gkk_intweight 204 elph_tr_ds%gamma_qpt_trout = elph_tr_ds%gamma_qpt_trout* elph_ds%occ_factor / elph_ds%k_phon%nkpt 205 elph_tr_ds%gamma_qpt_trin = elph_tr_ds%gamma_qpt_trin * elph_ds%occ_factor / elph_ds%k_phon%nkpt 206 207 write (message,'(2a)')' integrate_gamma_tr : transport gamma matrices are calculated ',& 208 & ' in recip space and for irred qpoints' 209 call wrtout(std_out,message,'COLL') 210 211 !DEBUG 212 !write(std_out,*)' integrate_gamma_tr_lova: end elph_tr_ds%gamma_qpt_trin(1,9,1,1,1)=',elph_tr_ds%gamma_qpt_trin(1,9,1,1,1) 213 !ENDDEBUG 214 215 216 end subroutine integrate_gamma_tr_lova