TABLE OF CONTENTS
ABINIT/nmsq_gam [ Functions ]
NAME
nmsq_gam
FUNCTION
Calculate gamma matrices keeping full dependence on bands from original h1_mat_el_sq matrix elements (no averaging over bands near the Fermi surface)
COPYRIGHT
Copyright (C) 2004-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
displ_red = phonon mode displacement vectors, post-multiplied by gprim matrix (ie. turned to reduced coordinates) eigvec = phonon eigenvectors elph_ds = datastructure with gkk matrix elements FSfullpqtofull = mapping of k+q to k kpt_phon = coordinates of kpoints near to FS h1_mat_el_sq = matrix elements $<psi_{k+q,m} | H^{1} | psi_{k,n}>$ squared iqptirred = index of present qpoint
OUTPUT
accum_mat = matrix for accumulating FS average of gkk (gamma matrix -> linewidths) accum_mat2 = matrix for accumulating FS average of gamma matrix with good prefactors
PARENTS
normsq_gkq
CHILDREN
gam_mult_displ,zgemm
SOURCE
41 #if defined HAVE_CONFIG_H 42 #include "config.h" 43 #endif 44 45 #include "abi_common.h" 46 47 48 subroutine nmsq_gam (accum_mat,accum_mat2,displ_red,eigvec,elph_ds,FSfullpqtofull,& 49 & h1_mat_el_sq,iqptirred) 50 51 use defs_basis 52 use defs_elphon 53 use m_errors 54 use m_profiling_abi 55 56 !This section has been created automatically by the script Abilint (TD). 57 !Do not modify the following lines by hand. 58 #undef ABI_FUNC 59 #define ABI_FUNC 'nmsq_gam' 60 use interfaces_77_ddb, except_this_one => nmsq_gam 61 !End of the abilint section 62 63 implicit none 64 65 !Arguments ------------------------------------ 66 !scalars 67 integer,intent(in) :: iqptirred 68 type(elph_type),intent(inout) :: elph_ds 69 !arrays 70 integer,intent(in) :: FSfullpqtofull(elph_ds%k_phon%nkpt,elph_ds%nqpt_full) 71 real(dp),intent(in) :: displ_red(2,elph_ds%nbranch,elph_ds%nbranch) 72 real(dp),intent(in) :: eigvec(2,elph_ds%nbranch,elph_ds%nbranch) 73 real(dp),intent(inout) :: & 74 & h1_mat_el_sq(2,elph_ds%nFSband*elph_ds%nFSband,elph_ds%nbranch*elph_ds%nbranch,elph_ds%k_phon%my_nkpt,elph_ds%nsppol) 75 real(dp),intent(inout) :: accum_mat(2,elph_ds%nbranch,elph_ds%nbranch,elph_ds%nsppol) 76 real(dp),intent(inout) :: accum_mat2(2,elph_ds%nbranch,elph_ds%nbranch,elph_ds%nsppol) 77 78 !Local variables------------------------------- 79 ! tmp variables for diagonalization 80 !scalars 81 integer :: ikpt_phon,ikpt_phonq,ib1,ib2,ibeff,ibranch,isppol,ipert1 82 integer :: jbranch 83 integer :: iqpt_fullbz 84 integer :: ik_this_proc 85 real(dp) :: sd1,sd2 86 character(len=500) :: message 87 !arrays 88 real(dp) :: gkq_1band(2,elph_ds%nbranch,elph_ds%nbranch) 89 real(dp) :: tmp_mat2(2,elph_ds%nbranch,elph_ds%nbranch) 90 real(dp) :: zgemm_tmp_mat(2,elph_ds%nbranch,elph_ds%nbranch) 91 92 ! ************************************************************************* 93 94 if (elph_ds%ep_keepbands == 0) then 95 write (message,'(a,i0)')' elph_ds%ep_keepbands should be 1 while is ',elph_ds%ep_keepbands 96 MSG_ERROR(message) 97 end if 98 99 !MG20060603 NOTE: 100 !accum_mat and accum_mat2 are real, the imaginary part is used for debugging purpose 101 !accum_mat2 is used to store the phonon-linewidhts before interpolation 102 103 iqpt_fullbz = elph_ds%qirredtofull(iqptirred) 104 write(std_out,*) 'nmsq_gam : iqptirred = ', iqptirred 105 106 do isppol=1,elph_ds%nsppol 107 do ik_this_proc =1, elph_ds%k_phon%my_nkpt 108 ikpt_phon = elph_ds%k_phon%my_ikpt(ik_this_proc) 109 110 ikpt_phonq = FSfullpqtofull(ikpt_phon,iqpt_fullbz) 111 112 do ib1=1,elph_ds%nFSband 113 sd1 = elph_ds%k_phon%wtk(ib1,ikpt_phon,isppol) !weights for distance from the fermi surface 114 115 do ib2=1,elph_ds%nFSband 116 sd2 = elph_ds%k_phon%wtk(ib2,ikpt_phonq,isppol) !weights for distance from the fermi surface 117 ibeff = ib2+elph_ds%nFSband*(ib1-1) 118 119 gkq_1band(:,:,:) = zero 120 121 zgemm_tmp_mat= reshape (h1_mat_el_sq(:,ibeff,:,ik_this_proc,isppol),(/2,elph_ds%nbranch,elph_ds%nbranch/)) 122 123 call gam_mult_displ(elph_ds%nbranch, displ_red, zgemm_tmp_mat, tmp_mat2) 124 125 ! sum over bands 126 do ipert1=1,elph_ds%nbranch 127 gkq_1band(1,ipert1,ipert1) = gkq_1band(1,ipert1,ipert1) + tmp_mat2(1,ipert1,ipert1) 128 end do 129 130 ! summing over k points and bands, still diagonal in jbranch 131 accum_mat(:,:,:,isppol) = accum_mat(:,:,:,isppol) + gkq_1band(:,:,:)*sd1*sd2 132 133 ! MG20060603 : summing over bands and kpoints with weights to calculate the phonon linewidth 134 do jbranch=1,elph_ds%nbranch 135 accum_mat2(:,jbranch,jbranch,isppol) = accum_mat2(:,jbranch,jbranch,isppol) + gkq_1band(:,jbranch,jbranch)*sd1*sd2 136 end do 137 ! END MG 138 139 140 ! now turn to cartesian coordinates 141 142 ! Final Gamma matrix (hermitian) = E * D_g * E^{+} 143 ! Where E^{+} is the hermitian conjugate of the eigenvector matrix E 144 ! And D_g is the diagonal matrix of values of gamma for this qpoint 145 146 ! Here gkq_1band is indexed with real phonon modes (not atom+idir) 147 ! turn gkq_1band to atom+cartesian coordinates (instead of normal coordinates for qpoint) 148 tmp_mat2(:,:,:) = zero 149 do ibranch =1,elph_ds%nbranch 150 do jbranch =1,elph_ds%nbranch 151 tmp_mat2(1,ibranch,jbranch) = tmp_mat2(1,ibranch,jbranch) + & 152 & eigvec(1,ibranch,jbranch) * gkq_1band(1,jbranch,jbranch) 153 tmp_mat2(2,ibranch,jbranch) = tmp_mat2(2,ibranch,jbranch) + & 154 & eigvec(2,ibranch,jbranch) * gkq_1band(1,jbranch,jbranch) 155 end do 156 end do 157 gkq_1band(:,:,:) = zero 158 159 ! here eigvec is transposed and complexconjugated. 160 zgemm_tmp_mat=zero 161 call zgemm('n','c',elph_ds%nbranch,elph_ds%nbranch,elph_ds%nbranch,cone,& 162 & tmp_mat2,elph_ds%nbranch,eigvec,elph_ds%nbranch,czero,zgemm_tmp_mat,elph_ds%nbranch) 163 164 gkq_1band = zgemm_tmp_mat 165 166 ! gamma matrix contribution in cartesian coordinates (ie interpolatable form) 167 h1_mat_el_sq(:,ibeff,:,ik_this_proc,isppol) = reshape(gkq_1band,(/2,elph_ds%nbranch*elph_ds%nbranch/)) 168 169 end do 170 end do 171 ! END loop over bands ib1 ib2 172 173 end do 174 ! END loop over kpt_phon 175 end do 176 !END loop over nsppol 177 178 179 end subroutine nmsq_gam