TABLE OF CONTENTS
ABINIT/prt_gkk_yambo [ Functions ]
NAME
prt_gkk_yambo
FUNCTION
This routine outputs el-phon related quantities for the yambo code at 1 q-point
COPYRIGHT
Copyright (C) 2009-2018 ABINIT group (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
displ_cart = phonon displacement vectors for this q-point in Cartesian coordinates. displ_red = phonon displacement vectors for this q-point, in reduced coordinates elph_ds = datastructure containing elphon matrix elements h1_mat_el = matrix elements of first order hamiltonian for present q-point, all perturbations iqptfull = index of present q-point in full array of q-points irredpert = index of irreducible perturbation (atom displaced) natom = number of atoms phfrq = phonon frequencies at present q-point qptn = q-point we will print for
OUTPUT
only writes to a file
NOTES
PARENTS
read_gkk
CHILDREN
SOURCE
43 #if defined HAVE_CONFIG_H 44 #include "config.h" 45 #endif 46 47 #include "abi_common.h" 48 49 50 subroutine prt_gkk_yambo(displ_cart,displ_red,kpt_phon,h1_mat_el,iqpt,& 51 & natom,nFSband,nkpt_phon,phfrq,qptn) 52 53 use defs_basis 54 use m_profiling_abi 55 56 use m_io_tools, only : get_unit 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 'prt_gkk_yambo' 62 !End of the abilint section 63 64 implicit none 65 66 !Arguments ------------------------------------ 67 !scalars 68 integer,intent(in) :: natom,iqpt 69 integer,intent(in) :: nFSband,nkpt_phon 70 !arrays 71 real(dp),intent(in) :: kpt_phon(3,nkpt_phon) 72 real(dp),intent(in) :: h1_mat_el(2,nFSband*nFSband,3*natom,nkpt_phon,1) 73 real(dp),intent(in) :: phfrq(3*natom) 74 real(dp),intent(in) :: displ_cart(2,3*natom,3*natom) 75 real(dp),intent(in) :: displ_red(2,3*natom,3*natom) 76 real(dp),intent(in) :: qptn(3) 77 78 !Local variables------------------------------- 79 !scalars 80 integer, save :: firsttime=1 81 integer :: outunit,ikpt,imode,iband,ibandp,iatom,idir,ibandindex 82 integer :: jmode, outunit2, outunit3 83 !arrays 84 real(dp) :: gkk_mode_dep(2) 85 ! ************************************************************************* 86 87 88 89 !if first time round: 90 if (firsttime==1) then 91 firsttime=0 92 ! squash file 93 outunit=get_unit() 94 open (unit=outunit,file="yambo_elphon_data",status="REPLACE") 95 outunit2=get_unit() 96 open (unit=outunit2,file="yambo_elphon_gkk_bymode",status="replace") 97 outunit3=get_unit() 98 open (unit=outunit3,file="yambo_elphon_gkksqtw_bymode",status="replace") 99 100 ! write dimensions 101 write (outunit,'(a,I6)') 'number of el atoms ', natom 102 write (outunit2,'(a,I6)') 'number of el atoms ', natom 103 write (outunit3,'(a,I6)') 'number of el atoms ', natom 104 write (outunit,'(a,I6)') 'number of ph modes ', 3*natom 105 write (outunit2,'(a,I6)') 'number of ph modes ', 3*natom 106 write (outunit3,'(a,I6)') 'number of ph modes ', 3*natom 107 write (outunit,'(a,I6)') 'number of el bands ', nFSband 108 write (outunit2,'(a,I6)') 'number of el bands ', nFSband 109 write (outunit3,'(a,I6)') 'number of el bands ', nFSband 110 111 ! write k-points 112 write (outunit,'(a,I6)') 'number of k-points ', nkpt_phon 113 write (outunit2,'(a,I6)') 'number of k-points ', nkpt_phon 114 write (outunit3,'(a,I6)') 'number of k-points ', nkpt_phon 115 do ikpt=1,nkpt_phon 116 write (outunit,'(a,I6,3E20.10)') 'reduced coord kpoint no ', ikpt, kpt_phon(:,ikpt) 117 write (outunit2,'(a,I6,3E20.10)') 'reduced coord kpoint no ', ikpt, kpt_phon(:,ikpt) 118 write (outunit3,'(a,I6,3E20.10)') 'reduced coord kpoint no ', ikpt, kpt_phon(:,ikpt) 119 end do 120 121 ! band energies are not accessible this deep in the code: simpler to get them 122 ! from elsewhere 123 124 close (outunit) 125 close (outunit2) 126 close (outunit3) 127 end if ! first time round 128 129 !open file 130 outunit=get_unit() 131 open (unit=outunit,file="yambo_elphon_data",status="unknown",position="append") 132 133 !qpoint 134 write (outunit,'(a,I6,3E20.10)') 'reduced coord qpoint no ', iqpt, qptn(:) 135 136 !frequencies 137 do imode=1,3*natom 138 write (outunit,'(a,I6,3E20.10)') 'phonon freq no ', imode, phfrq(imode) 139 end do 140 141 !displacement vector 142 do imode=1,3*natom 143 write (outunit,'(a,I6,3E20.10)') 'phonon displ vec no ', imode 144 do iatom=1,natom 145 write (outunit,'(3(2E20.10,2x))') displ_cart(:,(iatom-1)*3+1:iatom*3,imode) 146 end do 147 end do 148 149 !the beef: matrix elements of the first order hamiltonian for displacement of 150 !all atoms along all reduced directions 151 write (outunit,'(a)') ' matrix elements of all perturbations for this q-point' 152 do ikpt=1,nkpt_phon 153 write (outunit,'(a,I6)') ' kpoint ', ikpt 154 imode=0 155 do iatom=1,natom 156 do idir=1,3 157 imode=imode+1 158 write (outunit,'(a,I6,I6)') ' atom, direction = ', iatom,idir 159 ibandindex=0 160 do iband=1,nFSband 161 do ibandp=1,nFSband 162 ibandindex=ibandindex+1 163 write (outunit,'(a,I6,I6,2E20.10)') ' mat el for n,np ', iband,ibandp,& 164 & h1_mat_el(:,ibandindex,imode,ikpt,1) 165 end do !bandp 166 end do !band 167 end do !dir 168 end do !atom 169 end do 170 171 !blank line 172 write (outunit,*) 173 close (outunit) 174 175 outunit2=get_unit() 176 open (unit=outunit2,file="yambo_elphon_gkk_bymode",status="unknown",position="append") 177 outunit3=get_unit() 178 open (unit=outunit3,file="yambo_elphon_gkksqtw_bymode",status="unknown",position="append") 179 180 !qpoint 181 write (outunit2,'(a,I6,3E20.10)') 'reduced coord qpoint no ', iqpt, qptn(:) 182 write (outunit3,'(a,I6,3E20.10)') 'reduced coord qpoint no ', iqpt, qptn(:) 183 184 !print out mode-dependent matrix elements 185 write (outunit2,'(a)') ' matrix elements of all phonon modes for this q-point' 186 write (outunit3,'(a)') ' 1/w**1/2 times matrix elements of all phonon modes for this q-point' 187 do ikpt=1,nkpt_phon 188 write (outunit2,'(a,I6)') ' kpoint ', ikpt 189 write (outunit3,'(a,I6)') ' kpoint ', ikpt 190 ibandindex=0 191 do iband=1,nFSband 192 do ibandp=1,nFSband 193 ibandindex=ibandindex+1 194 write (outunit2,'(a,I6,I6)') ' el bands n,np ', iband,ibandp 195 write (outunit3,'(a,I6,I6)') ' el bands n,np ', iband,ibandp 196 do imode=1,3*natom 197 ! gkk_mode_dep = cg_zdotc(3*natom,displ_red(:,:,imode),h1_mat_el(:,ibandindex,:,ikpt,1)) 198 gkk_mode_dep = zero 199 do jmode=1,3*natom 200 gkk_mode_dep(1) = gkk_mode_dep(1) & 201 & + displ_red(1,jmode,imode)*h1_mat_el(1,ibandindex,jmode,ikpt,1) & 202 & + displ_red(2,jmode,imode)*h1_mat_el(2,ibandindex,jmode,ikpt,1) 203 gkk_mode_dep(2) = gkk_mode_dep(2) & 204 & + displ_red(1,jmode,imode)*h1_mat_el(2,ibandindex,jmode,ikpt,1) & 205 & - displ_red(2,jmode,imode)*h1_mat_el(1,ibandindex,jmode,ikpt,1) 206 end do 207 write (outunit2,'(a,I6,2E20.10)') ' mat el for phonon mode num = ', imode, gkk_mode_dep 208 write (outunit3,'(a,I6,2E20.10)') ' 1/w**1/2 * mat el for phonon mode num = ', & 209 & imode, gkk_mode_dep/sqrt(two*abs(phfrq(imode))+tol10) 210 end do !imode 211 end do !bandp 212 end do !band 213 end do 214 !blank line 215 write (outunit2,*) 216 write (outunit3,*) 217 218 close (outunit2) 219 close (outunit3) 220 221 end subroutine prt_gkk_yambo