TABLE OF CONTENTS
ABINIT/complete_gkk [ Functions ]
NAME
complete_gkk
FUNCTION
Use the set of special q points calculated by the Monkhorst & Pack Technique. Check if all the informations for the q points are present in the DDB to determine the elphon interaction matrices Generate the gkk matrices of the set of q points which samples homogeneously the entire Brillouin zone.
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
elph_ds = datastructure for elphon information (mainly matrix elements and dimensions) elph_ds%k_phon%full2full = kpt_phon index mapping under symops gkk_flag = flag for existence of matrix element gprimd(3,3)=dimensionful primitive translations in reciprocal space indsym = map of atoms by inverses of symrels natom=number of atoms in unit cell nsym=number of space group symmetries qpttoqpt = qpoint index mapping under symops rprimd(3,3)=dimensionful primitive translations in real space symrec(3,3,nsym)=3x3 matrices of the group symmetries (recip space) symrel(3,3,nsym)=3x3 matrices of the group symmetries (real space) tnons(3,nsym)=nonsymmorphic translations associated to symrel
OUTPUT
elph_ds%gkk_qpt = gkk matrices for all qpts on a full mesh
PARENTS
get_all_gkq
CHILDREN
xmpi_sum,zgemm
SOURCE
48 #if defined HAVE_CONFIG_H 49 #include "config.h" 50 #endif 51 52 #include "abi_common.h" 53 54 55 subroutine complete_gkk(elph_ds,gkk_flag,gprimd,indsym,natom,nsym,qpttoqpt,rprimd,symrec,symrel) 56 57 use defs_basis 58 use defs_elphon 59 use m_profiling_abi 60 use m_xmpi 61 62 !This section has been created automatically by the script Abilint (TD). 63 !Do not modify the following lines by hand. 64 #undef ABI_FUNC 65 #define ABI_FUNC 'complete_gkk' 66 !End of the abilint section 67 68 implicit none 69 70 !Arguments ------------------------------------ 71 !scalars 72 integer,intent(in) :: natom,nsym 73 type(elph_type),intent(inout) :: elph_ds 74 !arrays 75 integer,intent(in) :: indsym(4,nsym,natom) 76 integer,intent(in) :: qpttoqpt(2,nsym,elph_ds%nqpt_full),symrec(3,3,nsym) 77 integer,intent(in) :: symrel(3,3,nsym) 78 integer,intent(inout) :: gkk_flag(elph_ds%nbranch,elph_ds%nbranch,elph_ds%k_phon%my_nkpt,elph_ds%nsppol,elph_ds%nqpt_full) 79 real(dp),intent(in) :: gprimd(3,3) 80 real(dp),intent(in) :: rprimd(3,3) 81 82 !Local variables------------------------------- 83 !scalars 84 integer :: ikpt_phon,ib1,ibranch,ieqqpt,ii, ierr,comm 85 integer :: iqpt,isppol,isym 86 integer :: itim,jbranch,jj,kk,ll 87 integer :: neqqpt,symikpt_phon 88 integer :: iatom,ancestor_iatom 89 integer :: ik_this_proc, me,sz1,sz2 90 91 real(dp),parameter :: tol=2.d-8 92 !arrays 93 integer :: symmetrized_qpt(elph_ds%nqpt_full) 94 real(dp) :: ss(3,3) 95 real(dp) :: tmp_mat(2,elph_ds%nbranch,elph_ds%nbranch) 96 real(dp) :: tmp_mat2(2,elph_ds%nbranch,elph_ds%nbranch) 97 real(dp),allocatable :: gkk_qpt_new(:,:,:,:,:),gkk_qpt_tmp(:,:,:,:,:) 98 99 real(dp) :: ss_allatoms(2,elph_ds%nbranch,elph_ds%nbranch) 100 real(dp) :: c_one(2), c_zero(2) 101 102 103 ! ********************************************************************* 104 105 c_one = (/one,zero/) 106 c_zero = (/zero,zero/) 107 108 !Generation of the gkk matrices relative to the q points 109 !of the set which samples the entire Brillouin zone 110 111 comm = xmpi_world 112 me = xmpi_comm_rank(comm) 113 114 symmetrized_qpt(:) = -1 115 116 !FIXME bxu, why set it to 1? 117 !isppol=1 118 119 sz1=elph_ds%ngkkband*elph_ds%ngkkband 120 sz2=elph_ds%nbranch*elph_ds%nbranch 121 122 !these arrays are not parallelized, to enable symmetrization: syms swap k-points. 123 ABI_ALLOCATE(gkk_qpt_new,(2,sz1,sz2,elph_ds%k_phon%nkpt,elph_ds%nsppol)) 124 ABI_ALLOCATE(gkk_qpt_tmp,(2,sz1,sz2,elph_ds%k_phon%nkpt,elph_ds%nsppol)) 125 126 do iqpt=1,elph_ds%nqpt_full 127 128 ! Already symmetrized? 129 if (symmetrized_qpt(iqpt) == 1) cycle 130 131 gkk_qpt_new(:,:,:,:,:) = zero 132 ! gkk_qpt_tmp(:,:,:,:,:) = zero 133 134 ! loop over qpoints equivalent to iqpt 135 neqqpt=0 136 ! do not use time reversal symmetry to complete the qpoints: 137 ! do not know what happens to the gamma matrices 138 ! itim=1 139 140 do itim=1,2 141 do isym=1,nsym 142 ! ieqqpt is sent onto iqpt by itim/isym 143 ieqqpt = qpttoqpt(itim,isym,iqpt) 144 gkk_qpt_tmp(:,:,:,:,:) = zero 145 146 147 if (gkk_flag(1,1,1,1,ieqqpt) == -1) cycle 148 ! if we have information on this qpt 149 ! iqpt is equivalent to ieqqpt: get it from file or memory 150 do ik_this_proc =1,elph_ds%k_phon%my_nkpt 151 ikpt_phon = elph_ds%k_phon%my_ikpt(ik_this_proc) 152 153 if (elph_ds%gkqwrite == 0) then 154 gkk_qpt_tmp(:,:,:,ikpt_phon,:) = elph_ds%gkk_qpt(:,:,:,ik_this_proc,:,ieqqpt) 155 else if (elph_ds%gkqwrite == 1) then 156 read(elph_ds%unitgkq,REC=((ieqqpt-1)*elph_ds%k_phon%my_nkpt+ik_this_proc)) gkk_qpt_tmp(:,:,:,ikpt_phon,:) 157 end if 158 end do 159 160 ! condense everything 161 call xmpi_sum (gkk_qpt_tmp, comm, ierr) 162 163 neqqpt=neqqpt+1 164 165 if (elph_ds%ep_scalprod==1) then 166 do ii=1,3 167 do jj=1,3 168 ss(ii,jj)=0.0_dp 169 do kk=1,3 170 do ll=1,3 171 ss(ii,jj)=ss(ii,jj)+rprimd(ii,kk)*symrel(kk,ll,isym)*gprimd(ll,jj) 172 end do 173 end do 174 end do 175 end do 176 else 177 do ii=1,3 178 do jj=1,3 179 ss(ii,jj) = symrec(jj,ii,isym) 180 end do 181 end do 182 end if 183 184 ss_allatoms(:,:,:) = zero 185 do iatom=1,natom 186 ancestor_iatom = indsym(4,isym,iatom) 187 ! do jatom=1,natom 188 ! ancestor_jatom = indsym(4,isym,jatom) 189 ss_allatoms(1,(ancestor_iatom-1)*3+1:(ancestor_iatom-1)*3+3,& 190 & (iatom-1)*3+1: (iatom-1)*3+3) = ss(1:3,1:3) 191 ! end do 192 end do 193 194 195 ! NOTE ssinv(ii,jj)=ssinv(ii,jj)+gprimd(ii,kk)*rprimd(jj,ll)*symrec(ll,kk,isym) 196 197 do isppol=1,elph_ds%nsppol 198 do ikpt_phon=1,elph_ds%k_phon%nkpt 199 ! symikpt_phon is sent onto ikpt_phon by itim/isym 200 symikpt_phon=elph_ds%k_phon%full2full(itim,isym,ikpt_phon) 201 202 ! Do each element band1, band2 separately... 203 do ib1=1,elph_ds%ngkkband*elph_ds%ngkkband 204 205 ! multiply by the ss matrices 206 tmp_mat2(:,:,:) = zero 207 tmp_mat(:,:,:) = reshape(gkk_qpt_tmp(:,ib1,:,ikpt_phon,isppol),& 208 & (/2,elph_ds%nbranch,elph_ds%nbranch/)) 209 call ZGEMM ('N','N',elph_ds%nbranch,elph_ds%nbranch,elph_ds%nbranch,& 210 & c_one,ss_allatoms,elph_ds%nbranch,tmp_mat,elph_ds%nbranch,c_zero,& 211 & tmp_mat2,elph_ds%nbranch) 212 call ZGEMM ('N','T',elph_ds%nbranch,elph_ds%nbranch,elph_ds%nbranch,& 213 & c_one,tmp_mat2,elph_ds%nbranch,ss_allatoms,elph_ds%nbranch,c_zero,& 214 & tmp_mat,elph_ds%nbranch) 215 216 ! add to gkk_qpt_new 217 do ibranch =1,elph_ds%nbranch 218 do jbranch =1,elph_ds%nbranch 219 gkk_qpt_new(:,ib1,(jbranch-1)*elph_ds%nbranch+ibranch,symikpt_phon,isppol) = & 220 & gkk_qpt_new(:,ib1,(jbranch-1)*elph_ds%nbranch+ibranch,symikpt_phon,isppol) + & 221 & tmp_mat(:,jbranch,ibranch) 222 end do 223 end do 224 225 end do ! end ib1 do 226 end do ! end ikpt_phon do 227 end do ! end isppol do 228 229 end do ! end isym do 230 end do ! itim 231 232 if (neqqpt > 1) then 233 write(std_out,*) ' found several equiv qpts and am symmetrizing them ', neqqpt 234 end if 235 236 ! divide by number of equivalent qpts found 237 gkk_qpt_new(:,:,:,:,:) = gkk_qpt_new(:,:,:,:,:)/neqqpt 238 239 ! copy the symmetrized version into all the equivalent qpoints, appropriately transformed 240 ! See above 241 ! itim=1 242 do itim=1,2 243 do isym=1,nsym 244 ! ieqqpt is sent onto iqpt by itim/isym 245 ieqqpt = qpttoqpt(itim,isym,iqpt) 246 247 if (symmetrized_qpt(ieqqpt) /= -1) cycle 248 gkk_qpt_tmp(:,:,:,:,:) = zero 249 250 ! use symrec matrices to get inverse transform from isym^{-1} 251 if (elph_ds%ep_scalprod==1) then 252 do ii=1,3 253 do jj=1,3 254 ss(ii,jj)=0.0_dp 255 do kk=1,3 256 do ll=1,3 257 ! Use inverse of symop matrix here to get back to ieqqpt (inv+transpose is in symrec and in gprimd) 258 ss(ii,jj)=ss(ii,jj)+rprimd(ii,kk)*symrec(ll,kk,isym)*gprimd(ll,jj) 259 end do 260 end do 261 end do 262 end do 263 else 264 do ii=1,3 265 do jj=1,3 266 ss(ii,jj) = symrel(ii,jj,isym) 267 end do 268 end do 269 end if 270 271 ss_allatoms(:,:,:) = zero 272 do iatom=1,natom 273 ancestor_iatom = indsym(4,isym,iatom) 274 ! do jatom=1,natom 275 ! ancestor_jatom = indsym(4,isym,jatom) 276 ss_allatoms(1,(ancestor_iatom-1)*3+1:(ancestor_iatom-1)*3+3,& 277 & (iatom-1)*3+1: (iatom-1)*3+3) = ss(1:3,1:3) 278 ! end do 279 end do 280 281 ! ! Use inverse of symop matrix here to get back to ieqqpt 282 ! ssinv(ii,jj)=ssinv(ii,jj)+gprimd(ii,kk)*rprimd(jj,ll)*symrel(kk,ll,isym) 283 284 do isppol=1,elph_ds%nsppol 285 do ikpt_phon=1,elph_ds%k_phon%nkpt 286 ! symikpt_phon is sent onto ikpt_phon by itim/isym 287 symikpt_phon=elph_ds%k_phon%full2full(itim,isym,ikpt_phon) 288 289 do ib1=1,elph_ds%ngkkband*elph_ds%ngkkband 290 291 ! multiply by the ss^{-1} matrices 292 tmp_mat2(:,:,:) = zero 293 tmp_mat(:,:,:) = reshape(gkk_qpt_new(:,ib1,:,ikpt_phon,isppol),& 294 & (/2,elph_ds%nbranch,elph_ds%nbranch/)) 295 call ZGEMM ('N','N',elph_ds%nbranch,elph_ds%nbranch,elph_ds%nbranch,& 296 & c_one,ss_allatoms,elph_ds%nbranch,tmp_mat,elph_ds%nbranch,c_zero,& 297 & tmp_mat2,elph_ds%nbranch) 298 call ZGEMM ('N','T',elph_ds%nbranch,elph_ds%nbranch,elph_ds%nbranch,& 299 & c_one,tmp_mat2,elph_ds%nbranch,ss_allatoms,elph_ds%nbranch,c_zero,& 300 & tmp_mat,elph_ds%nbranch) 301 302 do ibranch =1,elph_ds%nbranch 303 do jbranch =1,elph_ds%nbranch 304 gkk_qpt_tmp(:,ib1,(jbranch-1)*elph_ds%nbranch+ibranch,symikpt_phon,isppol) =& 305 & tmp_mat(:,jbranch,ibranch) 306 end do 307 end do 308 309 do ik_this_proc =1,elph_ds%k_phon%my_nkpt 310 if (elph_ds%k_phon%my_ikpt(ik_this_proc) == symikpt_phon) then 311 if (gkk_flag (1,1,ik_this_proc,isppol,ieqqpt) == -1) gkk_flag (:,:,ik_this_proc,isppol,ieqqpt) = 0 312 exit 313 end if 314 end do 315 ! if (gkk_flag (1,1,symikpt_phon,isppol,ieqqpt) == -1) then 316 ! gkk_flag (:,:,symikpt_phon,isppol,ieqqpt) = 0 317 ! end if 318 319 end do ! end ib1 do 320 end do ! end ikpt_phon do 321 end do ! end isppol do 322 323 324 ! save symmetrized matrices for qpt ieqqpt 325 do ik_this_proc =1,elph_ds%k_phon%my_nkpt 326 ikpt_phon = elph_ds%k_phon%my_ikpt(ik_this_proc) 327 328 if (elph_ds%gkqwrite == 0) then 329 elph_ds%gkk_qpt(:,:,:,ik_this_proc,:,ieqqpt) = gkk_qpt_tmp(:,:,:,ikpt_phon,:) 330 else if (elph_ds%gkqwrite == 1) then 331 write(elph_ds%unitgkq,REC=((ieqqpt-1)*elph_ds%k_phon%my_nkpt+ik_this_proc)) gkk_qpt_tmp(:,:,:,ikpt_phon,:) 332 end if 333 end do 334 335 symmetrized_qpt(ieqqpt) = 1 336 337 end do ! end isym do 338 end do ! end itim do 339 340 end do 341 !end iqpt do 342 343 ABI_DEALLOCATE(gkk_qpt_new) 344 ABI_DEALLOCATE(gkk_qpt_tmp) 345 346 end subroutine complete_gkk