TABLE OF CONTENTS
ABINIT/read_el_veloc [ Functions ]
NAME
read_el_veloc
FUNCTION
This routine reads the velocities of the electronic GS for all kpts and bands then maps them into the FS kpt states
COPYRIGHT
Copyright (C) 2002-2018 ABINIT group (JPCroc) based on conducti This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt . For the initials of contributors, see ~abinit/doc/developers/contributors.txt .
INPUTS
nkpt_in = number of kpoints according to parent routine nband_in = number of bands according to parent routine nsppol_in = number of spin polarizations
OUTPUT
el_veloc(nkpt_in,nband_in,3)
PARENTS
elphon
CHILDREN
destroy_kptrank,get_rank_1kpt,hdr_free,inpgkk,mkkptrank
SOURCE
35 #if defined HAVE_CONFIG_H 36 #include "config.h" 37 #endif 38 39 #include "abi_common.h" 40 41 42 subroutine read_el_veloc(nband_in,nkpt_in,kpt_in,nsppol_in,elph_tr_ds) 43 44 use defs_basis 45 use defs_abitypes 46 use defs_elphon 47 use m_errors 48 use m_profiling_abi 49 use m_kptrank 50 use m_hdr 51 52 use m_io_tools, only : open_file 53 54 !This section has been created automatically by the script Abilint (TD). 55 !Do not modify the following lines by hand. 56 #undef ABI_FUNC 57 #define ABI_FUNC 'read_el_veloc' 58 use interfaces_72_response 59 !End of the abilint section 60 61 implicit none 62 63 !Arguments ----------------------------------- 64 !scalars 65 integer, intent(in) :: nband_in,nkpt_in,nsppol_in 66 type(elph_tr_type), intent(inout) :: elph_tr_ds 67 real(dp), intent(in) :: kpt_in(3,nkpt_in) 68 69 !Local variables------------------------------- 70 !scalars 71 integer :: bd2tot_index 72 integer :: iband,ii,ikpt, ikpt_ddk 73 integer :: isppol,l1,mband 74 integer :: bantot1 75 integer :: unit_ddk 76 integer :: symrankkpt 77 character(len=fnlen) :: filnam1,filnam2,filnam3 78 character(len=500) :: msg 79 type(hdr_type) :: hdr1 80 type(kptrank_type) :: kptrank_t 81 82 !arrays 83 real(dp) :: im_el_veloc(3) 84 real(dp),allocatable :: eig1_k(:,:) 85 real(dp),allocatable :: eigen11(:),eigen12(:),eigen13(:) 86 87 ! ********************************************************************************* 88 89 !Read data file name 90 !TODO: this should be standardized and read in anaddb always, not 91 !conditionally. Otherwise when new files are added to the anaddb files 92 !file... Catastrophe! 93 94 write(std_out,*)'enter read_el_veloc ' 95 96 !Read data file 97 if (open_file(elph_tr_ds%ddkfilename,msg,newunit=unit_ddk,form='formatted') /= 0) then 98 MSG_ERROR(msg) 99 end if 100 101 rewind(unit_ddk) 102 read(unit_ddk,'(a)')filnam1 ! first ddk file 103 read(unit_ddk,'(a)')filnam2 ! second ddk file 104 read(unit_ddk,'(a)')filnam3 ! third ddk file 105 close (unit_ddk) 106 107 bantot1 = 2*nband_in**2*nkpt_in*nsppol_in 108 109 call inpgkk(eigen11,filnam1,hdr1) 110 call hdr_free(hdr1) 111 112 call inpgkk(eigen12,filnam2,hdr1) 113 call hdr_free(hdr1) 114 115 !we use the hdr1 from the last call - should add some consistency 116 !testing here, we are trusting users not to mix different ddk files... 117 call inpgkk(eigen13,filnam3,hdr1) 118 119 !Extract info from the header 120 if(hdr1%nsppol /= nsppol_in) then 121 MSG_ERROR('nsspol /= input nsppol') 122 end if 123 124 !Get mband, as the maximum value of nband(nkpt) 125 mband=maxval(hdr1%nband(1:hdr1%nkpt)) 126 if (mband /= nband_in) then 127 MSG_ERROR('nband_in input to read_el_veloc is inconsistent with mband') 128 end if 129 130 write(std_out,*) 131 write(std_out,*) 'readings from read_el_veloc header' 132 write(std_out,'(a,i8)') ' natom =',hdr1%natom 133 write(std_out,'(a,3i8)') ' nkpt,nband_in,mband =',hdr1%nkpt,nband_in,mband 134 write(std_out,'(a, f10.5,a)' ) ' ecut =',hdr1%ecut,' Ha' 135 write(std_out,'(a,e15.5,a,e15.5,a)' )' fermie =',hdr1%fermie,' Ha ',hdr1%fermie*Ha_eV,' eV' 136 137 ABI_ALLOCATE(eig1_k,(2*nband_in**2,3)) 138 bd2tot_index = 0 139 elph_tr_ds%el_veloc=zero 140 141 !need correspondence between the DDK kpoints and the kpt_phon 142 call mkkptrank (hdr1%kptns,hdr1%nkpt,kptrank_t) 143 144 do isppol=1,nsppol_in 145 im_el_veloc(:)=zero 146 do ikpt=1,nkpt_in 147 call get_rank_1kpt (kpt_in(:,ikpt),symrankkpt, kptrank_t) 148 ikpt_ddk = kptrank_t%invrank(symrankkpt) 149 if (ikpt_ddk == -1) then 150 write(std_out,*)'read_el_veloc ******** error in correspondence between ddk and gkk kpoint sets' 151 write(std_out,*)' kpt sets in gkk and ddk files must agree.' 152 MSG_ERROR("Aborting now") 153 end if 154 bd2tot_index=2*nband_in**2*(ikpt_ddk-1) 155 156 ! first derivative eigenvalues for k-point 157 eig1_k(:,1)=eigen11(1+bd2tot_index:2*nband_in**2+bd2tot_index) 158 eig1_k(:,2)=eigen12(1+bd2tot_index:2*nband_in**2+bd2tot_index) 159 eig1_k(:,3)=eigen13(1+bd2tot_index:2*nband_in**2+bd2tot_index) 160 161 ! turn el_veloc to cartesian coordinates 162 do iband=1,nband_in 163 do l1=1,3 164 do ii=1,3 165 elph_tr_ds%el_veloc(ikpt,iband,l1,isppol)=elph_tr_ds%el_veloc(ikpt,iband,l1,isppol)+& 166 & hdr1%rprimd(l1,ii)*eig1_k(2*iband-1+(iband-1)*2*nband_in,ii)/two_pi 167 im_el_veloc(l1)=im_el_veloc(l1)+& 168 & hdr1%rprimd(l1,ii)*eig1_k(2*iband+(iband-1)*2*nband_in,ii)/two_pi 169 end do 170 end do ! l1 171 end do 172 end do 173 end do ! end isppol 174 175 call destroy_kptrank (kptrank_t) 176 ABI_DEALLOCATE(eig1_k) 177 ABI_DEALLOCATE(eigen11) 178 ABI_DEALLOCATE(eigen12) 179 ABI_DEALLOCATE(eigen13) 180 181 call hdr_free(hdr1) 182 183 write(std_out,*)'out of read_el_veloc ' 184 185 end subroutine read_el_veloc