TABLE OF CONTENTS
ABINIT/normsq_gkq [ Functions ]
NAME
normsq_gkq
FUNCTION
This routine takes the gkq matrix elements for a given qpoint, does the scalar product with the phonon displacement vector, squares the gkq matrix elements multiplies by the appropriate weights and puts them in a uniform (atom,icart) basis
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 in reduced coordinated. eigvec = eigenvectors of phonons (to turn to cartesian coord frame) elph_ds = datastructure with gkk matrix elements FSfullpqtofull = mapping of k + q to k h1_mat_el_sq = matrix elements $<psi_{k+q,m} | H^{1} | psi_{k,n}>$ matrix-squared iqptirred = index of present qpoint phfrq_tmp = phonon frequencies qpt_irred = array of qpoint coordinates
OUTPUT
elph_ds%gkq filled qdata(elph_ds%nbranch,elph_ds%nsppol,3) = array containing the phonon frequency, the linewidth and $\lambda_{q,\nu}$.
PARENTS
read_gkk
CHILDREN
gam_mult_displ,nmsq_gam,nmsq_gam_sumfs,nmsq_pure_gkk nmsq_pure_gkk_sumfs,wrtout,xmpi_sum,zhpev
SOURCE
43 #if defined HAVE_CONFIG_H 44 #include "config.h" 45 #endif 46 47 #include "abi_common.h" 48 49 subroutine normsq_gkq(displ_red,eigvec,elph_ds,FSfullpqtofull,& 50 & h1_mat_el_sq,iqptirred,phfrq_tmp,qpt_irred,qdata) 51 52 use defs_basis 53 use defs_elphon 54 use m_profiling_abi 55 use m_errors 56 use m_xmpi 57 58 !This section has been created automatically by the script Abilint (TD). 59 !Do not modify the following lines by hand. 60 #undef ABI_FUNC 61 #define ABI_FUNC 'normsq_gkq' 62 use interfaces_14_hidewrite 63 use interfaces_77_ddb, except_this_one => normsq_gkq 64 !End of the abilint section 65 66 implicit none 67 68 !Arguments ------------------------------------ 69 !scalars 70 integer,intent(in) :: iqptirred 71 type(elph_type),intent(inout) :: elph_ds 72 !arrays 73 integer,intent(in) :: FSfullpqtofull(elph_ds%k_phon%nkpt,elph_ds%nqpt_full) 74 real(dp),intent(in) :: displ_red(2,elph_ds%nbranch,elph_ds%nbranch) 75 real(dp),intent(in) :: eigvec(2,elph_ds%nbranch,elph_ds%nbranch) 76 real(dp),intent(inout) :: & 77 & 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) 78 real(dp),intent(in) :: phfrq_tmp(elph_ds%nbranch),qpt_irred(3,elph_ds%nqptirred) 79 real(dp),intent(out) :: qdata(elph_ds%nbranch,elph_ds%nsppol,3) 80 81 !Local variables------------------------------- 82 !scalars 83 integer :: i1,i2,ier,ii,isppol,jbranch,comm 84 real(dp) :: lambda_tot 85 character(len=500) :: message 86 !arrays 87 real(dp) :: accum_mat(2,elph_ds%nbranch,elph_ds%nbranch,elph_ds%nsppol) 88 real(dp) :: accum_mat2(2,elph_ds%nbranch,elph_ds%nbranch,elph_ds%nsppol) 89 real(dp) :: gam_now2(2,elph_ds%nbranch,elph_ds%nbranch) 90 real(dp) :: lambda(elph_ds%nsppol) 91 real(dp),allocatable :: matrx(:,:),val(:),vec(:,:,:) 92 real(dp),allocatable :: zhpev1(:,:),zhpev2(:) 93 94 ! ************************************************************************* 95 96 DBG_ENTER("COLL") 97 98 accum_mat = zero 99 accum_mat2 = zero 100 comm = xmpi_world 101 102 if (elph_ds%ep_scalprod == 1) then 103 ! 104 if (elph_ds%ep_keepbands == 0) then 105 call wrtout(std_out,' normsq_gkq : calling nmsq_gam_sumFS',"COLL") 106 call nmsq_gam_sumFS (accum_mat,accum_mat2,displ_red,eigvec,elph_ds,FSfullpqtofull,& 107 & h1_mat_el_sq,iqptirred) 108 109 else if (elph_ds%ep_keepbands == 1) then 110 call wrtout(std_out,' normsq_gkq : calling nmsq_gam',"COLL") 111 call nmsq_gam (accum_mat,accum_mat2,displ_red,eigvec,elph_ds,FSfullpqtofull,& 112 & h1_mat_el_sq,iqptirred) 113 114 else 115 write (message,'(a,i0)')' Wrong value for elph_ds%ep_keepbands = ',elph_ds%ep_keepbands 116 MSG_BUG(message) 117 end if 118 ! 119 else if (elph_ds%ep_scalprod == 0) then ! Interpolate on the pure "matrix of matrix elements" and do the scalar products later. 120 ! 121 if (elph_ds%ep_keepbands == 0) then 122 call wrtout(std_out,' normsq_gkq : calling nmsq_pure_gkk_sumFS',"COLL") 123 call nmsq_pure_gkk_sumFS (accum_mat,accum_mat2,displ_red,elph_ds,FSfullpqtofull,& 124 & h1_mat_el_sq,iqptirred) 125 126 else if (elph_ds%ep_keepbands == 1) then 127 call wrtout(std_out,' normsq_gkq : calling nmsq_pure_gkk',"COLL") 128 129 call nmsq_pure_gkk (accum_mat,accum_mat2,displ_red,elph_ds,FSfullpqtofull,& 130 & h1_mat_el_sq,iqptirred) 131 else 132 write (message,'(a,i0)')' Wrong value for elph_ds%ep_keepbands = ',elph_ds%ep_keepbands 133 MSG_BUG(message) 134 end if 135 ! 136 else 137 write (message,'(a,i0)')' Wrong value for elph_ds%ep_scalprod = ',elph_ds%ep_scalprod 138 MSG_BUG(message) 139 end if 140 !end if flag for doing scalar product now. 141 142 143 !MG: values without the good prefactor 144 accum_mat = accum_mat * elph_ds%occ_factor/elph_ds%k_phon%nkpt 145 146 !MG: accum_mat2 contains the line-widhts before the Fourier interpolation 147 accum_mat2 = accum_mat2 * elph_ds%occ_factor/elph_ds%k_phon%nkpt 148 149 !mpi sum over procs for accum_mat2 150 call xmpi_sum (accum_mat, comm, ier) 151 call xmpi_sum (accum_mat2, comm, ier) 152 153 !MG20060531i 154 !write e-ph quantities before Fourier interpolation 155 !save e-ph values in the temporary array qdata that will be copied into elph_ds%qgrid_data 156 157 write (message,'(4a,3es16.6,63a)')ch10, & 158 & ' Phonon linewidths before interpolation ',ch10, & 159 & ' Q point = ',qpt_irred(:,iqptirred),ch10,('=',ii=1,60),ch10,& 160 & ' Mode Frequency (Ha) Linewidth (Ha) Lambda ' 161 call wrtout(std_out,message,'COLL') 162 163 lambda_tot = zero 164 do isppol=1,elph_ds%nsppol 165 do ii=1,elph_ds%nbranch 166 lambda(isppol)=zero 167 ! MG: the tolerance factor is somehow arbitrary 168 if (abs(phfrq_tmp(ii)) > tol10) lambda(isppol)=accum_mat2(1,ii,ii,isppol)/& 169 & (pi*elph_ds%n0(isppol)*phfrq_tmp(ii)**2) 170 lambda_tot=lambda_tot+lambda(isppol) 171 write(message,'(i8,es20.6,2es16.6)' )ii,phfrq_tmp(ii),accum_mat2(1,ii,ii,isppol),lambda(isppol) 172 call wrtout(std_out,message,'COLL') 173 ! save values 174 qdata(ii,isppol,1)=phfrq_tmp(ii) 175 qdata(ii,isppol,2)=accum_mat2(1,ii,ii,isppol) 176 qdata(ii,isppol,3)=lambda(isppol) 177 end do !loop over branch 178 end do !loop over sppol 179 180 !normalize for number of spins 181 lambda_tot = lambda_tot / elph_ds%nsppol 182 183 write(message,'(61a,44x,es16.6,62a)' )('=',ii=1,60),ch10,lambda_tot,ch10,('=',ii=1,60),ch10 184 call wrtout(std_out,message,'COLL') 185 !ENDMG20060531 186 187 !immediately calculate linewidths: 188 write(std_out,*) 'summed accum_mat = ' 189 write(std_out,'(3(2E18.6,1x))') accum_mat(:,:,:,1) 190 write(std_out,*) 'summed accum_mat2 = ' 191 write(std_out,'(3(2E18.6,1x))') (accum_mat2(:,ii,ii,1),ii=1,elph_ds%nbranch) 192 write(std_out,*) 'displ_red = ' 193 write(std_out,'(3(2E18.6,1x))') displ_red 194 195 if (elph_ds%ep_scalprod == 1) then 196 do isppol=1,elph_ds%nsppol 197 ! Diagonalize gamma matrix at qpoint (complex matrix). Copied from dfpt_phfrq 198 ier=0 199 ii=1 200 ABI_ALLOCATE(matrx,(2,(elph_ds%nbranch*(elph_ds%nbranch+1))/2)) 201 do i2=1,elph_ds%nbranch 202 do i1=1,i2 203 matrx(1,ii)=accum_mat2(1,i1,i2,isppol) 204 matrx(2,ii)=accum_mat2(2,i1,i2,isppol) 205 ii=ii+1 206 end do 207 end do 208 ABI_ALLOCATE(zhpev1,(2,2*elph_ds%nbranch-1)) 209 ABI_ALLOCATE(zhpev2,(3*elph_ds%nbranch-2)) 210 ABI_ALLOCATE(val,(elph_ds%nbranch)) 211 ABI_ALLOCATE(vec,(2,elph_ds%nbranch,elph_ds%nbranch)) 212 call ZHPEV ('V','U',elph_ds%nbranch,matrx,val,vec,elph_ds%nbranch,zhpev1,zhpev2,ier) 213 214 write (std_out,*) ' normsq_gkq : accumulated eigenvalues isppol ',isppol, ' = ' 215 write (std_out,'(3E18.6)') val 216 ABI_DEALLOCATE(matrx) 217 ABI_DEALLOCATE(zhpev1) 218 ABI_DEALLOCATE(zhpev2) 219 ABI_DEALLOCATE(vec) 220 ABI_DEALLOCATE(val) 221 end do ! isppol 222 223 else if (elph_ds%ep_scalprod == 0) then 224 225 226 do isppol=1,elph_ds%nsppol 227 call gam_mult_displ(elph_ds%nbranch, displ_red, accum_mat(:,:,:,isppol), gam_now2) 228 229 write (std_out,*) ' normsq_gkq : accumulated eigenvalues isppol ', isppol, ' = ' 230 write (std_out,'(3(E14.6,1x))') (gam_now2(1,jbranch,jbranch), jbranch=1,elph_ds%nbranch) 231 write (std_out,*) ' normsq_gkq : imag part = ' 232 write (std_out,'(3(E14.6,1x))') (gam_now2(2,jbranch,jbranch), jbranch=1,elph_ds%nbranch) 233 end do ! isppol 234 235 end if 236 237 DBG_EXIT("COLL") 238 239 end subroutine normsq_gkq