TABLE OF CONTENTS
ABINIT/d2c_weights [ Functions ]
NAME
d2c_weights
FUNCTION
This routine calculates the integration weights on a coarse k-grid using the integration weights from a denser k-gird. The weights of the extra k points that being shared are evenly distributed. It also condenses the velocity*wtk and velcity^2*wtk.
COPYRIGHT
Copyright (C) 2012-2018 ABINIT group (BXU) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
elph_ds%k_fine%nkpt = number of fine FS k-points elph_ds%k_fine%wtk = integration weights of the fine FS k-grid elph_ds%k_phon%nkpt = number of coarse FS k-points elph_tr_ds%el_veloc = electronic velocities from the fine k-grid
OUTPUT
elph_ds%k_phon%wtk = integration weights of the coarse FS k-grid elph_ds%k_phon%velocwtk = velocity time integration weights of the coarse FS k-grid elph_ds%k_phon%vvelocwtk = velocity^2 time integration weights of the coarse FS k-grid
PARENTS
elphon,get_nv_fs_en
CHILDREN
SOURCE
36 #if defined HAVE_CONFIG_H 37 #include "config.h" 38 #endif 39 40 #include "abi_common.h" 41 42 subroutine d2c_weights(elph_ds,elph_tr_ds) 43 44 use defs_basis 45 use defs_elphon 46 use m_profiling_abi 47 use m_errors 48 49 !This section has been created automatically by the script Abilint (TD). 50 !Do not modify the following lines by hand. 51 #undef ABI_FUNC 52 #define ABI_FUNC 'd2c_weights' 53 !End of the abilint section 54 55 implicit none 56 57 !Arguments ------------------------------------ 58 type(elph_type),intent(inout) :: elph_ds 59 type(elph_tr_type),intent(inout),optional :: elph_tr_ds 60 61 !Local variables------------------------------- 62 integer :: ii, jj, kk 63 integer :: ikpt, jkpt, kkpt 64 integer :: iikpt, jjkpt, kkkpt 65 integer :: icomp, jcomp 66 integer :: iFSband, iband 67 integer :: ikpt_fine, ikpt_phon 68 integer :: nkpt_fine1, nkpt_phon1 69 integer :: nkpt_fine2, nkpt_phon2 70 integer :: nkpt_fine3, nkpt_phon3 71 integer :: nscale1, nscale2, nscale3 72 73 ! ************************************************************************* 74 nkpt_phon1 = elph_ds%kptrlatt(1,1) 75 nkpt_phon2 = elph_ds%kptrlatt(2,2) 76 nkpt_phon3 = elph_ds%kptrlatt(3,3) 77 nkpt_fine1 = elph_ds%kptrlatt_fine(1,1) 78 nkpt_fine2 = elph_ds%kptrlatt_fine(2,2) 79 nkpt_fine3 = elph_ds%kptrlatt_fine(3,3) 80 nscale1 = dble(nkpt_fine1/nkpt_phon1) 81 nscale2 = dble(nkpt_fine2/nkpt_phon2) 82 nscale3 = dble(nkpt_fine3/nkpt_phon3) 83 if (abs(INT(nscale1)-nscale1) > 0.01) then 84 MSG_ERROR('The denser k-gird MUST be multiples of the phon k-grid') 85 end if 86 if (abs(INT(nscale2)-nscale2) > 0.01) then 87 MSG_ERROR('The denser k-gird MUST be multiples of the phon k-grid') 88 end if 89 if (abs(INT(nscale3)-nscale3) > 0.01) then 90 MSG_ERROR('The denser k-gird MUST be multiples of the phon k-grid') 91 end if 92 nscale1 = INT(nscale1) 93 nscale2 = INT(nscale2) 94 nscale3 = INT(nscale3) 95 96 !bxu, get wtk of coarse grid from fine grid 97 elph_ds%k_phon%wtk = zero 98 if (present(elph_tr_ds)) then 99 elph_ds%k_phon%velocwtk = zero 100 elph_ds%k_phon%vvelocwtk = zero 101 end if 102 103 do ikpt = 1, nkpt_phon1 104 do jkpt = 1, nkpt_phon2 105 do kkpt = 1, nkpt_phon3 106 ikpt_phon = kkpt + (jkpt-1)*nkpt_phon3 + (ikpt-1)*nkpt_phon2*nkpt_phon3 107 ! inside the paralellepipe 108 do ii = -((nscale1+1)/2-1), ((nscale1+1)/2-1) 109 do jj = -((nscale2+1)/2-1), ((nscale2+1)/2-1) 110 do kk = -((nscale3+1)/2-1), ((nscale3+1)/2-1) 111 iikpt = 1 + (ikpt-1)*nscale1 + ii 112 jjkpt = 1 + (jkpt-1)*nscale2 + jj 113 kkkpt = 1 + (kkpt-1)*nscale3 + kk 114 if (iikpt .le. 0) iikpt = iikpt + nkpt_fine1 115 if (jjkpt .le. 0) jjkpt = jjkpt + nkpt_fine2 116 if (kkkpt .le. 0) kkkpt = kkkpt + nkpt_fine3 117 if (iikpt .gt. nkpt_fine1) iikpt = iikpt - nkpt_fine1 118 if (jjkpt .gt. nkpt_fine2) jjkpt = jjkpt - nkpt_fine2 119 if (kkkpt .gt. nkpt_fine3) kkkpt = kkkpt - nkpt_fine3 120 ikpt_fine = kkkpt + (jjkpt-1)*nkpt_fine3 + (iikpt-1)*nkpt_fine2*nkpt_fine3 121 elph_ds%k_phon%wtk(:,ikpt_phon,:)=elph_ds%k_phon%wtk(:,ikpt_phon,:)+& 122 & elph_ds%k_fine%wtk(:,ikpt_fine,:) 123 if (present(elph_tr_ds)) then 124 do iFSband=1,elph_ds%ngkkband !FS bands 125 iband=iFSband+elph_ds%minFSband-1 ! full bands 126 do icomp = 1, 3 127 elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)=elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)+& 128 & elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:) 129 do jcomp = 1, 3 130 elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:) = & 131 & elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:)+& 132 & elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)* & 133 & elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,jcomp,:) 134 end do 135 end do 136 end do 137 end if 138 end do 139 end do 140 end do 141 ! on the 6 faces 142 if (MOD(nscale3,2) == 0) then ! when nscale3 is an even number 143 do ii = -((nscale1+1)/2-1), ((nscale1+1)/2-1) 144 do jj = -((nscale2+1)/2-1), ((nscale2+1)/2-1) 145 iikpt = 1 + (ikpt-1)*nscale1 + ii 146 jjkpt = 1 + (jkpt-1)*nscale2 + jj 147 if (iikpt .le. 0) iikpt = iikpt + nkpt_fine1 148 if (jjkpt .le. 0) jjkpt = jjkpt + nkpt_fine2 149 if (iikpt .gt. nkpt_fine1) iikpt = iikpt - nkpt_fine1 150 if (jjkpt .gt. nkpt_fine2) jjkpt = jjkpt - nkpt_fine2 151 152 kkkpt = 1 + (kkpt-1)*nscale3 + nscale3/2 153 if (kkkpt .le. 0) kkkpt = kkkpt + nkpt_fine3 154 if (kkkpt .gt. nkpt_fine3) kkkpt = kkkpt - nkpt_fine3 155 ikpt_fine = kkkpt + (jjkpt-1)*nkpt_fine3 + (iikpt-1)*nkpt_fine2*nkpt_fine3 156 elph_ds%k_phon%wtk(:,ikpt_phon,:)=elph_ds%k_phon%wtk(:,ikpt_phon,:)+& 157 & 0.5_dp*elph_ds%k_fine%wtk(:,ikpt_fine,:) 158 if (present(elph_tr_ds)) then 159 do iFSband=1,elph_ds%ngkkband !FS bands 160 iband=iFSband+elph_ds%minFSband-1 ! full bands 161 do icomp = 1, 3 162 elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)=elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)+& 163 & 0.5_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:) 164 do jcomp = 1, 3 165 elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:) = & 166 & elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:)+& 167 & 0.5_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)* & 168 & elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,jcomp,:) 169 end do 170 end do 171 end do 172 end if 173 174 kkkpt = 1 + (kkpt-1)*nscale3 - nscale3/2 175 if (kkkpt .le. 0) kkkpt = kkkpt + nkpt_fine3 176 if (kkkpt .gt. nkpt_fine3) kkkpt = kkkpt - nkpt_fine3 177 ikpt_fine = kkkpt + (jjkpt-1)*nkpt_fine3 + (iikpt-1)*nkpt_fine2*nkpt_fine3 178 elph_ds%k_phon%wtk(:,ikpt_phon,:)=elph_ds%k_phon%wtk(:,ikpt_phon,:)+& 179 & 0.5_dp*elph_ds%k_fine%wtk(:,ikpt_fine,:) 180 if (present(elph_tr_ds)) then 181 do iFSband=1,elph_ds%ngkkband !FS bands 182 iband=iFSband+elph_ds%minFSband-1 ! full bands 183 do icomp = 1, 3 184 elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)=elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)+& 185 & 0.5_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:) 186 do jcomp = 1, 3 187 elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:) = & 188 & elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:)+& 189 & 0.5_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)* & 190 & elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,jcomp,:) 191 end do 192 end do 193 end do 194 end if 195 end do 196 end do 197 end if 198 if (MOD(nscale2,2) == 0) then ! when nscale2 is an even number 199 do ii = -((nscale1+1)/2-1), ((nscale1+1)/2-1) 200 do kk = -((nscale3+1)/2-1), ((nscale3+1)/2-1) 201 iikpt = 1 + (ikpt-1)*nscale1 + ii 202 kkkpt = 1 + (kkpt-1)*nscale3 + kk 203 if (iikpt .le. 0) iikpt = iikpt + nkpt_fine1 204 if (kkkpt .le. 0) kkkpt = kkkpt + nkpt_fine3 205 if (iikpt .gt. nkpt_fine1) iikpt = iikpt - nkpt_fine1 206 if (kkkpt .gt. nkpt_fine3) kkkpt = kkkpt - nkpt_fine3 207 208 jjkpt = 1 + (jkpt-1)*nscale2 + nscale2/2 209 if (jjkpt .le. 0) jjkpt = jjkpt + nkpt_fine2 210 if (jjkpt .gt. nkpt_fine2) jjkpt = jjkpt - nkpt_fine2 211 ikpt_fine = kkkpt + (jjkpt-1)*nkpt_fine3 + (iikpt-1)*nkpt_fine2*nkpt_fine3 212 elph_ds%k_phon%wtk(:,ikpt_phon,:)=elph_ds%k_phon%wtk(:,ikpt_phon,:)+& 213 & 0.5_dp*elph_ds%k_fine%wtk(:,ikpt_fine,:) 214 if (present(elph_tr_ds)) then 215 do iFSband=1,elph_ds%ngkkband !FS bands 216 iband=iFSband+elph_ds%minFSband-1 ! full bands 217 do icomp = 1, 3 218 elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)=elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)+& 219 & 0.5_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:) 220 do jcomp = 1, 3 221 elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:) = & 222 & elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:)+& 223 & 0.5_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)* & 224 & elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,jcomp,:) 225 end do 226 end do 227 end do 228 end if 229 230 jjkpt = 1 + (jkpt-1)*nscale2 - nscale2/2 231 if (jjkpt .le. 0) jjkpt = jjkpt + nkpt_fine2 232 if (jjkpt .gt. nkpt_fine2) jjkpt = jjkpt - nkpt_fine2 233 ikpt_fine = kkkpt + (jjkpt-1)*nkpt_fine3 + (iikpt-1)*nkpt_fine2*nkpt_fine3 234 elph_ds%k_phon%wtk(:,ikpt_phon,:)=elph_ds%k_phon%wtk(:,ikpt_phon,:)+& 235 & 0.5_dp*elph_ds%k_fine%wtk(:,ikpt_fine,:) 236 if (present(elph_tr_ds)) then 237 do iFSband=1,elph_ds%ngkkband !FS bands 238 iband=iFSband+elph_ds%minFSband-1 ! full bands 239 do icomp = 1, 3 240 elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)=elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)+& 241 & 0.5_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:) 242 do jcomp = 1, 3 243 elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:) = & 244 & elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:)+& 245 & 0.5_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)* & 246 & elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,jcomp,:) 247 end do 248 end do 249 end do 250 end if 251 end do 252 end do 253 end if 254 if (MOD(nscale1,2) == 0) then ! when nscale1 is an even number 255 do kk = -((nscale3+1)/2-1), ((nscale3+1)/2-1) 256 do jj = -((nscale2+1)/2-1), ((nscale2+1)/2-1) 257 kkkpt = 1 + (kkpt-1)*nscale3 + kk 258 jjkpt = 1 + (jkpt-1)*nscale2 + jj 259 if (kkkpt .le. 0) kkkpt = kkkpt + nkpt_fine3 260 if (jjkpt .le. 0) jjkpt = jjkpt + nkpt_fine2 261 if (kkkpt .gt. nkpt_fine3) kkkpt = kkkpt - nkpt_fine3 262 if (jjkpt .gt. nkpt_fine2) jjkpt = jjkpt - nkpt_fine2 263 264 iikpt = 1 + (ikpt-1)*nscale1 + nscale1/2 265 if (iikpt .le. 0) iikpt = iikpt + nkpt_fine1 266 if (iikpt .gt. nkpt_fine1) iikpt = iikpt - nkpt_fine1 267 ikpt_fine = kkkpt + (jjkpt-1)*nkpt_fine3 + (iikpt-1)*nkpt_fine2*nkpt_fine3 268 elph_ds%k_phon%wtk(:,ikpt_phon,:)=elph_ds%k_phon%wtk(:,ikpt_phon,:)+& 269 & 0.5_dp*elph_ds%k_fine%wtk(:,ikpt_fine,:) 270 if (present(elph_tr_ds)) then 271 do iFSband=1,elph_ds%ngkkband !FS bands 272 iband=iFSband+elph_ds%minFSband-1 ! full bands 273 do icomp = 1, 3 274 elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)=elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)+& 275 & 0.5_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:) 276 do jcomp = 1, 3 277 elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:) = & 278 & elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:)+& 279 & 0.5_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)* & 280 & elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,jcomp,:) 281 end do 282 end do 283 end do 284 end if 285 286 iikpt = 1 + (ikpt-1)*nscale1 - nscale1/2 287 if (iikpt .le. 0) iikpt = iikpt + nkpt_fine1 288 if (iikpt .gt. nkpt_fine1) iikpt = iikpt - nkpt_fine1 289 ikpt_fine = kkkpt + (jjkpt-1)*nkpt_fine3 + (iikpt-1)*nkpt_fine2*nkpt_fine3 290 elph_ds%k_phon%wtk(:,ikpt_phon,:)=elph_ds%k_phon%wtk(:,ikpt_phon,:)+& 291 & 0.5_dp*elph_ds%k_fine%wtk(:,ikpt_fine,:) 292 if (present(elph_tr_ds)) then 293 do iFSband=1,elph_ds%ngkkband !FS bands 294 iband=iFSband+elph_ds%minFSband-1 ! full bands 295 do icomp = 1, 3 296 elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)=elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)+& 297 & 0.5_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:) 298 do jcomp = 1, 3 299 elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:) = & 300 & elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:)+& 301 & 0.5_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)* & 302 & elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,jcomp,:) 303 end do 304 end do 305 end do 306 end if 307 end do 308 end do 309 ! on the 12 sides 310 end if 311 if (MOD(nscale2,2) == 0 .and. MOD(nscale3,2) == 0) then 312 do ii = -((nscale1+1)/2-1), ((nscale1+1)/2-1) 313 iikpt = 1 + (ikpt-1)*nscale1 + ii 314 if (iikpt .le. 0) iikpt = iikpt + nkpt_fine1 315 if (iikpt .gt. nkpt_fine1) iikpt = iikpt - nkpt_fine1 316 317 jjkpt = 1 + (jkpt-1)*nscale2 + nscale2/2 318 kkkpt = 1 + (kkpt-1)*nscale3 + nscale3/2 319 if (jjkpt .le. 0) jjkpt = jjkpt + nkpt_fine2 320 if (jjkpt .gt. nkpt_fine2) jjkpt = jjkpt - nkpt_fine2 321 if (kkkpt .le. 0) kkkpt = kkkpt + nkpt_fine3 322 if (kkkpt .gt. nkpt_fine3) kkkpt = kkkpt - nkpt_fine3 323 ikpt_fine = kkkpt + (jjkpt-1)*nkpt_fine3 + (iikpt-1)*nkpt_fine2*nkpt_fine3 324 elph_ds%k_phon%wtk(:,ikpt_phon,:)=elph_ds%k_phon%wtk(:,ikpt_phon,:)+& 325 & 0.25_dp*elph_ds%k_fine%wtk(:,ikpt_fine,:) 326 if (present(elph_tr_ds)) then 327 do iFSband=1,elph_ds%ngkkband !FS bands 328 iband=iFSband+elph_ds%minFSband-1 ! full bands 329 do icomp = 1, 3 330 elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)=elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)+& 331 & 0.25_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:) 332 do jcomp = 1, 3 333 elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:) = & 334 & elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:)+& 335 & 0.25_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)* & 336 & elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,jcomp,:) 337 end do 338 end do 339 end do 340 end if 341 342 jjkpt = 1 + (jkpt-1)*nscale2 + nscale2/2 343 kkkpt = 1 + (kkpt-1)*nscale3 - nscale3/2 344 if (jjkpt .le. 0) jjkpt = jjkpt + nkpt_fine2 345 if (jjkpt .gt. nkpt_fine2) jjkpt = jjkpt - nkpt_fine2 346 if (kkkpt .le. 0) kkkpt = kkkpt + nkpt_fine3 347 if (kkkpt .gt. nkpt_fine3) kkkpt = kkkpt - nkpt_fine3 348 ikpt_fine = kkkpt + (jjkpt-1)*nkpt_fine3 + (iikpt-1)*nkpt_fine2*nkpt_fine3 349 elph_ds%k_phon%wtk(:,ikpt_phon,:)=elph_ds%k_phon%wtk(:,ikpt_phon,:)+& 350 & 0.25_dp*elph_ds%k_fine%wtk(:,ikpt_fine,:) 351 if (present(elph_tr_ds)) then 352 do iFSband=1,elph_ds%ngkkband !FS bands 353 iband=iFSband+elph_ds%minFSband-1 ! full bands 354 do icomp = 1, 3 355 elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)=elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)+& 356 & 0.25_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:) 357 do jcomp = 1, 3 358 elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:) = & 359 & elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:)+& 360 & 0.25_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)* & 361 & elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,jcomp,:) 362 end do 363 end do 364 end do 365 end if 366 367 jjkpt = 1 + (jkpt-1)*nscale2 - nscale2/2 368 kkkpt = 1 + (kkpt-1)*nscale3 + nscale3/2 369 if (jjkpt .le. 0) jjkpt = jjkpt + nkpt_fine2 370 if (jjkpt .gt. nkpt_fine2) jjkpt = jjkpt - nkpt_fine2 371 if (kkkpt .le. 0) kkkpt = kkkpt + nkpt_fine3 372 if (kkkpt .gt. nkpt_fine3) kkkpt = kkkpt - nkpt_fine3 373 ikpt_fine = kkkpt + (jjkpt-1)*nkpt_fine3 + (iikpt-1)*nkpt_fine2*nkpt_fine3 374 elph_ds%k_phon%wtk(:,ikpt_phon,:)=elph_ds%k_phon%wtk(:,ikpt_phon,:)+& 375 & 0.25_dp*elph_ds%k_fine%wtk(:,ikpt_fine,:) 376 if (present(elph_tr_ds)) then 377 do iFSband=1,elph_ds%ngkkband !FS bands 378 iband=iFSband+elph_ds%minFSband-1 ! full bands 379 do icomp = 1, 3 380 elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)=elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)+& 381 & 0.25_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:) 382 do jcomp = 1, 3 383 elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:) = & 384 & elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:)+& 385 & 0.25_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)* & 386 & elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,jcomp,:) 387 end do 388 end do 389 end do 390 end if 391 392 jjkpt = 1 + (jkpt-1)*nscale2 - nscale2/2 393 kkkpt = 1 + (kkpt-1)*nscale3 - nscale3/2 394 if (jjkpt .le. 0) jjkpt = jjkpt + nkpt_fine2 395 if (jjkpt .gt. nkpt_fine2) jjkpt = jjkpt - nkpt_fine2 396 if (kkkpt .le. 0) kkkpt = kkkpt + nkpt_fine3 397 if (kkkpt .gt. nkpt_fine3) kkkpt = kkkpt - nkpt_fine3 398 ikpt_fine = kkkpt + (jjkpt-1)*nkpt_fine3 + (iikpt-1)*nkpt_fine2*nkpt_fine3 399 elph_ds%k_phon%wtk(:,ikpt_phon,:)=elph_ds%k_phon%wtk(:,ikpt_phon,:)+& 400 & 0.25_dp*elph_ds%k_fine%wtk(:,ikpt_fine,:) 401 if (present(elph_tr_ds)) then 402 do iFSband=1,elph_ds%ngkkband !FS bands 403 iband=iFSband+elph_ds%minFSband-1 ! full bands 404 do icomp = 1, 3 405 elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)=elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)+& 406 & 0.25_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:) 407 do jcomp = 1, 3 408 elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:) = & 409 & elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:)+& 410 & 0.25_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)* & 411 & elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,jcomp,:) 412 end do 413 end do 414 end do 415 end if 416 end do 417 end if 418 if (MOD(nscale1,2) == 0 .and. MOD(nscale3,2) == 0) then 419 do jj = -((nscale2+1)/2-1), ((nscale2+1)/2-1) 420 jjkpt = 1 + (jkpt-1)*nscale2 + jj 421 if (jjkpt .le. 0) jjkpt = jjkpt + nkpt_fine2 422 if (jjkpt .gt. nkpt_fine2) jjkpt = jjkpt - nkpt_fine2 423 424 iikpt = 1 + (ikpt-1)*nscale1 + nscale1/2 425 kkkpt = 1 + (kkpt-1)*nscale3 + nscale3/2 426 if (iikpt .le. 0) iikpt = iikpt + nkpt_fine1 427 if (iikpt .gt. nkpt_fine1) iikpt = iikpt - nkpt_fine1 428 if (kkkpt .le. 0) kkkpt = kkkpt + nkpt_fine3 429 if (kkkpt .gt. nkpt_fine3) kkkpt = kkkpt - nkpt_fine3 430 ikpt_fine = kkkpt + (jjkpt-1)*nkpt_fine3 + (iikpt-1)*nkpt_fine2*nkpt_fine3 431 elph_ds%k_phon%wtk(:,ikpt_phon,:)=elph_ds%k_phon%wtk(:,ikpt_phon,:)+& 432 & 0.25_dp*elph_ds%k_fine%wtk(:,ikpt_fine,:) 433 if (present(elph_tr_ds)) then 434 do iFSband=1,elph_ds%ngkkband !FS bands 435 iband=iFSband+elph_ds%minFSband-1 ! full bands 436 do icomp = 1, 3 437 elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)=elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)+& 438 & 0.25_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:) 439 do jcomp = 1, 3 440 elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:) = & 441 & elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:)+& 442 & 0.25_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)* & 443 & elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,jcomp,:) 444 end do 445 end do 446 end do 447 end if 448 449 iikpt = 1 + (ikpt-1)*nscale1 + nscale1/2 450 kkkpt = 1 + (kkpt-1)*nscale3 - nscale3/2 451 if (iikpt .le. 0) iikpt = iikpt + nkpt_fine1 452 if (iikpt .gt. nkpt_fine1) iikpt = iikpt - nkpt_fine1 453 if (kkkpt .le. 0) kkkpt = kkkpt + nkpt_fine3 454 if (kkkpt .gt. nkpt_fine3) kkkpt = kkkpt - nkpt_fine3 455 ikpt_fine = kkkpt + (jjkpt-1)*nkpt_fine3 + (iikpt-1)*nkpt_fine2*nkpt_fine3 456 elph_ds%k_phon%wtk(:,ikpt_phon,:)=elph_ds%k_phon%wtk(:,ikpt_phon,:)+& 457 & 0.25_dp*elph_ds%k_fine%wtk(:,ikpt_fine,:) 458 if (present(elph_tr_ds)) then 459 do iFSband=1,elph_ds%ngkkband !FS bands 460 iband=iFSband+elph_ds%minFSband-1 ! full bands 461 do icomp = 1, 3 462 elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)=elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)+& 463 & 0.25_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:) 464 do jcomp = 1, 3 465 elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:) = & 466 & elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:)+& 467 & 0.25_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)* & 468 & elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,jcomp,:) 469 end do 470 end do 471 end do 472 end if 473 474 iikpt = 1 + (ikpt-1)*nscale1 - nscale1/2 475 kkkpt = 1 + (kkpt-1)*nscale3 + nscale3/2 476 if (iikpt .le. 0) iikpt = iikpt + nkpt_fine1 477 if (iikpt .gt. nkpt_fine1) iikpt = iikpt - nkpt_fine1 478 if (kkkpt .le. 0) kkkpt = kkkpt + nkpt_fine3 479 if (kkkpt .gt. nkpt_fine3) kkkpt = kkkpt - nkpt_fine3 480 ikpt_fine = kkkpt + (jjkpt-1)*nkpt_fine3 + (iikpt-1)*nkpt_fine2*nkpt_fine3 481 elph_ds%k_phon%wtk(:,ikpt_phon,:)=elph_ds%k_phon%wtk(:,ikpt_phon,:)+& 482 & 0.25_dp*elph_ds%k_fine%wtk(:,ikpt_fine,:) 483 if (present(elph_tr_ds)) then 484 do iFSband=1,elph_ds%ngkkband !FS bands 485 iband=iFSband+elph_ds%minFSband-1 ! full bands 486 do icomp = 1, 3 487 elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)=elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)+& 488 & 0.25_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:) 489 do jcomp = 1, 3 490 elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:) = & 491 & elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:)+& 492 & 0.25_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)* & 493 & elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,jcomp,:) 494 end do 495 end do 496 end do 497 end if 498 499 iikpt = 1 + (ikpt-1)*nscale1 - nscale1/2 500 kkkpt = 1 + (kkpt-1)*nscale3 - nscale3/2 501 if (iikpt .le. 0) iikpt = iikpt + nkpt_fine1 502 if (iikpt .gt. nkpt_fine1) iikpt = iikpt - nkpt_fine1 503 if (kkkpt .le. 0) kkkpt = kkkpt + nkpt_fine3 504 if (kkkpt .gt. nkpt_fine3) kkkpt = kkkpt - nkpt_fine3 505 ikpt_fine = kkkpt + (jjkpt-1)*nkpt_fine3 + (iikpt-1)*nkpt_fine2*nkpt_fine3 506 elph_ds%k_phon%wtk(:,ikpt_phon,:)=elph_ds%k_phon%wtk(:,ikpt_phon,:)+& 507 & 0.25_dp*elph_ds%k_fine%wtk(:,ikpt_fine,:) 508 if (present(elph_tr_ds)) then 509 do iFSband=1,elph_ds%ngkkband !FS bands 510 iband=iFSband+elph_ds%minFSband-1 ! full bands 511 do icomp = 1, 3 512 elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)=elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)+& 513 & 0.25_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:) 514 do jcomp = 1, 3 515 elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:) = & 516 & elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:)+& 517 & 0.25_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)* & 518 & elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,jcomp,:) 519 end do 520 end do 521 end do 522 end if 523 end do 524 end if 525 if (MOD(nscale2,2) == 0 .and. MOD(nscale1,2) == 0) then 526 do kk = -((nscale3+1)/2-1), ((nscale3+1)/2-1) 527 kkkpt = 1 + (kkpt-1)*nscale3 + kk 528 if (kkkpt .le. 0) kkkpt = kkkpt + nkpt_fine3 529 if (kkkpt .gt. nkpt_fine3) kkkpt = kkkpt - nkpt_fine3 530 531 jjkpt = 1 + (jkpt-1)*nscale2 + nscale2/2 532 iikpt = 1 + (ikpt-1)*nscale1 + nscale1/2 533 if (jjkpt .le. 0) jjkpt = jjkpt + nkpt_fine2 534 if (jjkpt .gt. nkpt_fine2) jjkpt = jjkpt - nkpt_fine2 535 if (iikpt .le. 0) iikpt = iikpt + nkpt_fine1 536 if (iikpt .gt. nkpt_fine1) iikpt = iikpt - nkpt_fine1 537 ikpt_fine = kkkpt + (jjkpt-1)*nkpt_fine3 + (iikpt-1)*nkpt_fine2*nkpt_fine3 538 elph_ds%k_phon%wtk(:,ikpt_phon,:)=elph_ds%k_phon%wtk(:,ikpt_phon,:)+& 539 & 0.25_dp*elph_ds%k_fine%wtk(:,ikpt_fine,:) 540 if (present(elph_tr_ds)) then 541 do iFSband=1,elph_ds%ngkkband !FS bands 542 iband=iFSband+elph_ds%minFSband-1 ! full bands 543 do icomp = 1, 3 544 elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)=elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)+& 545 & 0.25_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:) 546 do jcomp = 1, 3 547 elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:) = & 548 & elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:)+& 549 & 0.25_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)* & 550 & elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,jcomp,:) 551 end do 552 end do 553 end do 554 end if 555 556 jjkpt = 1 + (jkpt-1)*nscale2 + nscale2/2 557 iikpt = 1 + (ikpt-1)*nscale1 - nscale1/2 558 if (jjkpt .le. 0) jjkpt = jjkpt + nkpt_fine2 559 if (jjkpt .gt. nkpt_fine2) jjkpt = jjkpt - nkpt_fine2 560 if (iikpt .le. 0) iikpt = iikpt + nkpt_fine1 561 if (iikpt .gt. nkpt_fine1) iikpt = iikpt - nkpt_fine1 562 ikpt_fine = kkkpt + (jjkpt-1)*nkpt_fine3 + (iikpt-1)*nkpt_fine2*nkpt_fine3 563 elph_ds%k_phon%wtk(:,ikpt_phon,:)=elph_ds%k_phon%wtk(:,ikpt_phon,:)+& 564 & 0.25_dp*elph_ds%k_fine%wtk(:,ikpt_fine,:) 565 if (present(elph_tr_ds)) then 566 do iFSband=1,elph_ds%ngkkband !FS bands 567 iband=iFSband+elph_ds%minFSband-1 ! full bands 568 do icomp = 1, 3 569 elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)=elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)+& 570 & 0.25_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:) 571 do jcomp = 1, 3 572 elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:) = & 573 & elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:)+& 574 & 0.25_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)* & 575 & elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,jcomp,:) 576 end do 577 end do 578 end do 579 end if 580 581 jjkpt = 1 + (jkpt-1)*nscale2 - nscale2/2 582 iikpt = 1 + (ikpt-1)*nscale1 + nscale1/2 583 if (jjkpt .le. 0) jjkpt = jjkpt + nkpt_fine2 584 if (jjkpt .gt. nkpt_fine2) jjkpt = jjkpt - nkpt_fine2 585 if (iikpt .le. 0) iikpt = iikpt + nkpt_fine1 586 if (iikpt .gt. nkpt_fine1) iikpt = iikpt - nkpt_fine1 587 ikpt_fine = kkkpt + (jjkpt-1)*nkpt_fine3 + (iikpt-1)*nkpt_fine2*nkpt_fine3 588 elph_ds%k_phon%wtk(:,ikpt_phon,:)=elph_ds%k_phon%wtk(:,ikpt_phon,:)+& 589 & 0.25_dp*elph_ds%k_fine%wtk(:,ikpt_fine,:) 590 if (present(elph_tr_ds)) then 591 do iFSband=1,elph_ds%ngkkband !FS bands 592 iband=iFSband+elph_ds%minFSband-1 ! full bands 593 do icomp = 1, 3 594 elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)=elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)+& 595 & 0.25_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:) 596 do jcomp = 1, 3 597 elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:) = & 598 & elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:)+& 599 & 0.25_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)* & 600 & elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,jcomp,:) 601 end do 602 end do 603 end do 604 end if 605 606 jjkpt = 1 + (jkpt-1)*nscale2 - nscale2/2 607 iikpt = 1 + (ikpt-1)*nscale1 - nscale1/2 608 if (jjkpt .le. 0) jjkpt = jjkpt + nkpt_fine2 609 if (jjkpt .gt. nkpt_fine2) jjkpt = jjkpt - nkpt_fine2 610 if (iikpt .le. 0) iikpt = iikpt + nkpt_fine1 611 if (iikpt .gt. nkpt_fine1) iikpt = iikpt - nkpt_fine1 612 ikpt_fine = kkkpt + (jjkpt-1)*nkpt_fine3 + (iikpt-1)*nkpt_fine2*nkpt_fine3 613 elph_ds%k_phon%wtk(:,ikpt_phon,:)=elph_ds%k_phon%wtk(:,ikpt_phon,:)+& 614 & 0.25_dp*elph_ds%k_fine%wtk(:,ikpt_fine,:) 615 if (present(elph_tr_ds)) then 616 do iFSband=1,elph_ds%ngkkband !FS bands 617 iband=iFSband+elph_ds%minFSband-1 ! full bands 618 do icomp = 1, 3 619 elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)=elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)+& 620 & 0.25_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:) 621 do jcomp = 1, 3 622 elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:) = & 623 & elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:)+& 624 & 0.25_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)* & 625 & elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,jcomp,:) 626 end do 627 end do 628 end do 629 end if 630 end do 631 ! on the 8 corners 632 end if 633 if (MOD(nscale1,2) == 0 .and. MOD(nscale2,2) == 0 .and. MOD(nscale3,2) == 0) then 634 iikpt = 1 + (ikpt-1)*nscale1 + nscale1/2 635 jjkpt = 1 + (jkpt-1)*nscale2 + nscale2/2 636 kkkpt = 1 + (kkpt-1)*nscale3 + nscale3/2 637 if (iikpt .le. 0) iikpt = iikpt + nkpt_fine1 638 if (jjkpt .le. 0) jjkpt = jjkpt + nkpt_fine2 639 if (kkkpt .le. 0) kkkpt = kkkpt + nkpt_fine3 640 if (iikpt .gt. nkpt_fine1) iikpt = iikpt - nkpt_fine1 641 if (jjkpt .gt. nkpt_fine2) jjkpt = jjkpt - nkpt_fine2 642 if (kkkpt .gt. nkpt_fine3) kkkpt = kkkpt - nkpt_fine3 643 ikpt_fine = kkkpt + (jjkpt-1)*nkpt_fine3 + (iikpt-1)*nkpt_fine2*nkpt_fine3 644 elph_ds%k_phon%wtk(:,ikpt_phon,:)=elph_ds%k_phon%wtk(:,ikpt_phon,:)+& 645 & 0.125_dp*elph_ds%k_fine%wtk(:,ikpt_fine,:) 646 if (present(elph_tr_ds)) then 647 do iFSband=1,elph_ds%ngkkband !FS bands 648 iband=iFSband+elph_ds%minFSband-1 ! full bands 649 do icomp = 1, 3 650 elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)=elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)+& 651 & 0.125_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:) 652 do jcomp = 1, 3 653 elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:) = & 654 & elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:)+& 655 & 0.125_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)* & 656 & elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,jcomp,:) 657 end do 658 end do 659 end do 660 end if 661 662 iikpt = 1 + (ikpt-1)*nscale1 + nscale1/2 663 jjkpt = 1 + (jkpt-1)*nscale2 + nscale2/2 664 kkkpt = 1 + (kkpt-1)*nscale3 - nscale3/2 665 if (iikpt .le. 0) iikpt = iikpt + nkpt_fine1 666 if (jjkpt .le. 0) jjkpt = jjkpt + nkpt_fine2 667 if (kkkpt .le. 0) kkkpt = kkkpt + nkpt_fine3 668 if (iikpt .gt. nkpt_fine1) iikpt = iikpt - nkpt_fine1 669 if (jjkpt .gt. nkpt_fine2) jjkpt = jjkpt - nkpt_fine2 670 if (kkkpt .gt. nkpt_fine3) kkkpt = kkkpt - nkpt_fine3 671 ikpt_fine = kkkpt + (jjkpt-1)*nkpt_fine3 + (iikpt-1)*nkpt_fine2*nkpt_fine3 672 elph_ds%k_phon%wtk(:,ikpt_phon,:)=elph_ds%k_phon%wtk(:,ikpt_phon,:)+& 673 & 0.125_dp*elph_ds%k_fine%wtk(:,ikpt_fine,:) 674 if (present(elph_tr_ds)) then 675 do iFSband=1,elph_ds%ngkkband !FS bands 676 iband=iFSband+elph_ds%minFSband-1 ! full bands 677 do icomp = 1, 3 678 elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)=elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)+& 679 & 0.125_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:) 680 do jcomp = 1, 3 681 elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:) = & 682 & elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:)+& 683 & 0.125_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)* & 684 & elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,jcomp,:) 685 end do 686 end do 687 end do 688 end if 689 690 iikpt = 1 + (ikpt-1)*nscale1 + nscale1/2 691 jjkpt = 1 + (jkpt-1)*nscale2 - nscale2/2 692 kkkpt = 1 + (kkpt-1)*nscale3 + nscale3/2 693 if (iikpt .le. 0) iikpt = iikpt + nkpt_fine1 694 if (jjkpt .le. 0) jjkpt = jjkpt + nkpt_fine2 695 if (kkkpt .le. 0) kkkpt = kkkpt + nkpt_fine3 696 if (iikpt .gt. nkpt_fine1) iikpt = iikpt - nkpt_fine1 697 if (jjkpt .gt. nkpt_fine2) jjkpt = jjkpt - nkpt_fine2 698 if (kkkpt .gt. nkpt_fine3) kkkpt = kkkpt - nkpt_fine3 699 ikpt_fine = kkkpt + (jjkpt-1)*nkpt_fine3 + (iikpt-1)*nkpt_fine2*nkpt_fine3 700 elph_ds%k_phon%wtk(:,ikpt_phon,:)=elph_ds%k_phon%wtk(:,ikpt_phon,:)+& 701 & 0.125_dp*elph_ds%k_fine%wtk(:,ikpt_fine,:) 702 if (present(elph_tr_ds)) then 703 do iFSband=1,elph_ds%ngkkband !FS bands 704 iband=iFSband+elph_ds%minFSband-1 ! full bands 705 do icomp = 1, 3 706 elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)=elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)+& 707 & 0.125_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:) 708 do jcomp = 1, 3 709 elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:) = & 710 & elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:)+& 711 & 0.125_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)* & 712 & elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,jcomp,:) 713 end do 714 end do 715 end do 716 end if 717 718 iikpt = 1 + (ikpt-1)*nscale1 + nscale1/2 719 jjkpt = 1 + (jkpt-1)*nscale2 - nscale2/2 720 kkkpt = 1 + (kkpt-1)*nscale3 - nscale3/2 721 if (iikpt .le. 0) iikpt = iikpt + nkpt_fine1 722 if (jjkpt .le. 0) jjkpt = jjkpt + nkpt_fine2 723 if (kkkpt .le. 0) kkkpt = kkkpt + nkpt_fine3 724 if (iikpt .gt. nkpt_fine1) iikpt = iikpt - nkpt_fine1 725 if (jjkpt .gt. nkpt_fine2) jjkpt = jjkpt - nkpt_fine2 726 if (kkkpt .gt. nkpt_fine3) kkkpt = kkkpt - nkpt_fine3 727 ikpt_fine = kkkpt + (jjkpt-1)*nkpt_fine3 + (iikpt-1)*nkpt_fine2*nkpt_fine3 728 elph_ds%k_phon%wtk(:,ikpt_phon,:)=elph_ds%k_phon%wtk(:,ikpt_phon,:)+& 729 & 0.125_dp*elph_ds%k_fine%wtk(:,ikpt_fine,:) 730 if (present(elph_tr_ds)) then 731 do iFSband=1,elph_ds%ngkkband !FS bands 732 iband=iFSband+elph_ds%minFSband-1 ! full bands 733 do icomp = 1, 3 734 elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)=elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)+& 735 & 0.125_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:) 736 do jcomp = 1, 3 737 elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:) = & 738 & elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:)+& 739 & 0.125_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)* & 740 & elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,jcomp,:) 741 end do 742 end do 743 end do 744 end if 745 746 iikpt = 1 + (ikpt-1)*nscale1 - nscale1/2 747 jjkpt = 1 + (jkpt-1)*nscale2 + nscale2/2 748 kkkpt = 1 + (kkpt-1)*nscale3 + nscale3/2 749 if (iikpt .le. 0) iikpt = iikpt + nkpt_fine1 750 if (jjkpt .le. 0) jjkpt = jjkpt + nkpt_fine2 751 if (kkkpt .le. 0) kkkpt = kkkpt + nkpt_fine3 752 if (iikpt .gt. nkpt_fine1) iikpt = iikpt - nkpt_fine1 753 if (jjkpt .gt. nkpt_fine2) jjkpt = jjkpt - nkpt_fine2 754 if (kkkpt .gt. nkpt_fine3) kkkpt = kkkpt - nkpt_fine3 755 ikpt_fine = kkkpt + (jjkpt-1)*nkpt_fine3 + (iikpt-1)*nkpt_fine2*nkpt_fine3 756 elph_ds%k_phon%wtk(:,ikpt_phon,:)=elph_ds%k_phon%wtk(:,ikpt_phon,:)+& 757 & 0.125_dp*elph_ds%k_fine%wtk(:,ikpt_fine,:) 758 if (present(elph_tr_ds)) then 759 do iFSband=1,elph_ds%ngkkband !FS bands 760 iband=iFSband+elph_ds%minFSband-1 ! full bands 761 do icomp = 1, 3 762 elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)=elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)+& 763 & 0.125_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:) 764 do jcomp = 1, 3 765 elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:) = & 766 & elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:)+& 767 & 0.125_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)* & 768 & elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,jcomp,:) 769 end do 770 end do 771 end do 772 end if 773 774 iikpt = 1 + (ikpt-1)*nscale1 - nscale1/2 775 jjkpt = 1 + (jkpt-1)*nscale2 + nscale2/2 776 kkkpt = 1 + (kkpt-1)*nscale3 - nscale3/2 777 if (iikpt .le. 0) iikpt = iikpt + nkpt_fine1 778 if (jjkpt .le. 0) jjkpt = jjkpt + nkpt_fine2 779 if (kkkpt .le. 0) kkkpt = kkkpt + nkpt_fine3 780 if (iikpt .gt. nkpt_fine1) iikpt = iikpt - nkpt_fine1 781 if (jjkpt .gt. nkpt_fine2) jjkpt = jjkpt - nkpt_fine2 782 if (kkkpt .gt. nkpt_fine3) kkkpt = kkkpt - nkpt_fine3 783 ikpt_fine = kkkpt + (jjkpt-1)*nkpt_fine3 + (iikpt-1)*nkpt_fine2*nkpt_fine3 784 elph_ds%k_phon%wtk(:,ikpt_phon,:)=elph_ds%k_phon%wtk(:,ikpt_phon,:)+& 785 & 0.125_dp*elph_ds%k_fine%wtk(:,ikpt_fine,:) 786 if (present(elph_tr_ds)) then 787 do iFSband=1,elph_ds%ngkkband !FS bands 788 iband=iFSband+elph_ds%minFSband-1 ! full bands 789 do icomp = 1, 3 790 elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)=elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)+& 791 & 0.125_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:) 792 do jcomp = 1, 3 793 elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:) = & 794 & elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:)+& 795 & 0.125_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)* & 796 & elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,jcomp,:) 797 end do 798 end do 799 end do 800 end if 801 802 iikpt = 1 + (ikpt-1)*nscale1 - nscale1/2 803 jjkpt = 1 + (jkpt-1)*nscale2 - nscale2/2 804 kkkpt = 1 + (kkpt-1)*nscale3 + nscale3/2 805 if (iikpt .le. 0) iikpt = iikpt + nkpt_fine1 806 if (jjkpt .le. 0) jjkpt = jjkpt + nkpt_fine2 807 if (kkkpt .le. 0) kkkpt = kkkpt + nkpt_fine3 808 if (iikpt .gt. nkpt_fine1) iikpt = iikpt - nkpt_fine1 809 if (jjkpt .gt. nkpt_fine2) jjkpt = jjkpt - nkpt_fine2 810 if (kkkpt .gt. nkpt_fine3) kkkpt = kkkpt - nkpt_fine3 811 ikpt_fine = kkkpt + (jjkpt-1)*nkpt_fine3 + (iikpt-1)*nkpt_fine2*nkpt_fine3 812 elph_ds%k_phon%wtk(:,ikpt_phon,:)=elph_ds%k_phon%wtk(:,ikpt_phon,:)+& 813 & 0.125_dp*elph_ds%k_fine%wtk(:,ikpt_fine,:) 814 if (present(elph_tr_ds)) then 815 do iFSband=1,elph_ds%ngkkband !FS bands 816 iband=iFSband+elph_ds%minFSband-1 ! full bands 817 do icomp = 1, 3 818 elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)=elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)+& 819 & 0.125_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:) 820 do jcomp = 1, 3 821 elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:) = & 822 & elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:)+& 823 & 0.125_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)* & 824 & elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,jcomp,:) 825 end do 826 end do 827 end do 828 end if 829 830 iikpt = 1 + (ikpt-1)*nscale1 - nscale1/2 831 jjkpt = 1 + (jkpt-1)*nscale2 - nscale2/2 832 kkkpt = 1 + (kkpt-1)*nscale3 - nscale3/2 833 if (iikpt .le. 0) iikpt = iikpt + nkpt_fine1 834 if (jjkpt .le. 0) jjkpt = jjkpt + nkpt_fine2 835 if (kkkpt .le. 0) kkkpt = kkkpt + nkpt_fine3 836 if (iikpt .gt. nkpt_fine1) iikpt = iikpt - nkpt_fine1 837 if (jjkpt .gt. nkpt_fine2) jjkpt = jjkpt - nkpt_fine2 838 if (kkkpt .gt. nkpt_fine3) kkkpt = kkkpt - nkpt_fine3 839 ikpt_fine = kkkpt + (jjkpt-1)*nkpt_fine3 + (iikpt-1)*nkpt_fine2*nkpt_fine3 840 elph_ds%k_phon%wtk(:,ikpt_phon,:)=elph_ds%k_phon%wtk(:,ikpt_phon,:)+& 841 & 0.125_dp*elph_ds%k_fine%wtk(:,ikpt_fine,:) 842 if (present(elph_tr_ds)) then 843 do iFSband=1,elph_ds%ngkkband !FS bands 844 iband=iFSband+elph_ds%minFSband-1 ! full bands 845 do icomp = 1, 3 846 elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)=elph_ds%k_phon%velocwtk(iFSband,ikpt_phon,icomp,:)+& 847 & 0.125_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:) 848 do jcomp = 1, 3 849 elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:) = & 850 & elph_ds%k_phon%vvelocwtk(iFSband,ikpt_phon,icomp,jcomp,:)+& 851 & 0.125_dp*elph_ds%k_fine%wtk(iFSband,ikpt_fine,:)* & 852 & elph_tr_ds%el_veloc(ikpt_fine,iband,icomp,:)*elph_tr_ds%el_veloc(ikpt_fine,iband,jcomp,:) 853 end do 854 end do 855 end do 856 end if 857 end if 858 end do 859 end do 860 end do 861 862 !bxu, divide by nscale^3 to be consistent with the normalization of kpt_phon 863 elph_ds%k_phon%wtk = elph_ds%k_phon%wtk/nscale1/nscale2/nscale3 864 if (present(elph_tr_ds)) then 865 elph_ds%k_phon%velocwtk = elph_ds%k_phon%velocwtk/nscale1/nscale2/nscale3 866 elph_ds%k_phon%vvelocwtk = elph_ds%k_phon%vvelocwtk/nscale1/nscale2/nscale3 867 end if 868 869 end subroutine d2c_weights