TABLE OF CONTENTS
ABINIT/dfptff_gbefd [ Functions ]
NAME
dfptff_gbefd
FUNCTION
calculate the gradient of the second order \Omega E \cdot P term, Eq.(23) in PRB 75, 115116(2007).
COPYRIGHT
Copyright (C) 2004-2018 ABINIT group (XW). 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 .
INPUTS
cg(2,mpw*nspinor*mband*mkmem*nsppol) = planewave coefficients of wavefunctions cg1(2,mpw1*nspinor*mband*mk1mem*nsppol) = pw coefficients of RF wavefunctions at k,q. dtefield = variables related to response Berry-phase calculation ikpt = the index of the current k point isppol = the index of the spin component mband = maximum number of bands mkmem = maximum number of k-points in core memory mpw = maximum number of plane waves mpw1 = maximum number of plane waves for response wavefunctions nkpt = number of k points npwarr(nkpt) = number of planewaves in basis and boundary at this k point npwar1(nkpt) = number of planewaves in basis and boundary for response wfs nspinor = 1 for scalar wfs, 2 for spinor wfs nsppol = 1 for unpolarized, 2 for spin-polarized qmat(2,dtefield%mband_occ,dtefield%mband_occ,nkpt,2,3) = inverse of the overlap matrix pwindall(max(mpw,mpw1)*mkmem,8,3) = array used to compute the overlap matrices pwindall(:,1,:) <- <u^(0)_i|u^(0)_i+1> pwindall(:,2,:) <- <u^(0)_i|u^(0)_i-1> pwindall(:,3,:) <- <u^(1)_i|u^(1)_i+1> pwindall(:,4,:) <- <u^(1)_i|u^(1)_i-1> pwindall(:,5,:) <- <u^(1)_i|u^(0)_i+n+1> pwindall(:,6,:) <- <u^(1)_i|u^(0)_i+n-1> pwindall(:,7,:) <- <u^(0)_i|u^(1)_i-n+1> pwindall(:,8,:) <- <u^(0)_i|u^(1)_i-n-1>
OUTPUT
grad_berry = the gradient of the Berry phase term
PARENTS
dfpt_vtorho
CHILDREN
overlap_g
SOURCE
56 #if defined HAVE_CONFIG_H 57 #include "config.h" 58 #endif 59 60 #include "abi_common.h" 61 62 63 subroutine dfptff_gbefd(cg,cg1,dtefield,grad_berry,idir_efield,ikpt,isppol,mband,mpw,mpw1,mkmem,mk1mem,nkpt,& 64 & npwarr,npwar1,nspinor,nsppol,qmat,pwindall,rprimd) 65 66 67 use defs_basis 68 use m_profiling_abi 69 use m_efield 70 71 use m_cgtools, only : overlap_g 72 73 !This section has been created automatically by the script Abilint (TD). 74 !Do not modify the following lines by hand. 75 #undef ABI_FUNC 76 #define ABI_FUNC 'dfptff_gbefd' 77 !End of the abilint section 78 79 implicit none 80 81 !Arguments ---------------------------------------- 82 !scalars 83 integer,intent(in) :: idir_efield,ikpt,isppol,mband,mk1mem,mkmem,mpw,mpw1,nkpt 84 integer,intent(in) :: nspinor,nsppol 85 type(efield_type),intent(in) :: dtefield 86 !arrays 87 integer,intent(in) :: npwar1(nkpt),npwarr(nkpt) 88 integer,intent(in) :: pwindall(max(mpw,mpw1)*mkmem,8,3) 89 real(dp),intent(in) :: cg(2,mpw*nspinor*mband*mkmem*nsppol) 90 real(dp),intent(in) :: cg1(2,mpw1*nspinor*mband*mk1mem*nsppol) 91 real(dp),intent(in) :: qmat(2,dtefield%mband_occ,dtefield%mband_occ,nkpt,2,3) 92 real(dp),intent(in) :: rprimd(3,3) 93 real(dp),intent(out) :: grad_berry(2,mpw1,dtefield%mband_occ) 94 95 !Local variables ------------------------- 96 !scalars 97 integer :: iband,icg,icg1,idir,ikpt1 98 integer :: ikptn,ikptnp1,ipw,jband,jpw,kband 99 integer :: mpw_tmp,npw_k1,npw_k2,pwmax,pwmin 100 real(dp) :: doti,dotr,fac,wfi,wfr 101 !arrays 102 integer,allocatable :: pwind_tmp(:) 103 real(dp) :: z1(2),z2(2) 104 real(dp),allocatable :: Amat(:,:,:),Bmat(:,:,:),s1mat(:,:,:),vect1(:,:) 105 real(dp),allocatable :: vect2(:,:) 106 107 ! ************************************************************************* 108 109 mpw_tmp=max(mpw,mpw1) 110 ABI_ALLOCATE(vect1,(2,0:mpw_tmp)) 111 ABI_ALLOCATE(vect2,(2,0:mpw_tmp)) 112 ABI_ALLOCATE(s1mat,(2,dtefield%mband_occ,dtefield%mband_occ)) 113 ABI_ALLOCATE(pwind_tmp,(mpw_tmp)) 114 ABI_ALLOCATE(Amat,(2,dtefield%mband_occ,dtefield%mband_occ)) 115 ABI_ALLOCATE(Bmat,(2,dtefield%mband_occ,dtefield%mband_occ)) 116 vect1(:,0) = zero ; vect2(:,0) = zero 117 s1mat(:,:,:)=zero 118 grad_berry(:,:,:) = zero 119 120 do idir=1,3 121 fac = dtefield%efield_dot(idir)*dble(nkpt)/& 122 & (dble(dtefield%nstr(idir))*four_pi) 123 124 ! prepare 125 ikpt1 = dtefield%ikpt_dk(ikpt,1,idir) 126 icg1 = dtefield%cgindex(ikpt1,isppol+nsppol) 127 npw_k1 = npwar1(ikpt) 128 npw_k2 = npwar1(ikpt1) 129 pwind_tmp(1:npw_k1) = pwindall((ikpt-1)*mpw_tmp+1:(ikpt-1)*mpw_tmp+npw_k1,3,idir) 130 131 do ipw = 1, npw_k1 132 jpw = pwind_tmp(ipw) 133 if (jpw > 0) then 134 do iband = 1, dtefield%mband_occ 135 wfr = cg1(1,icg1 + (iband - 1)*npw_k2*nspinor + jpw) 136 wfi = cg1(2,icg1 + (iband - 1)*npw_k2*nspinor + jpw) 137 do jband = 1, dtefield%mband_occ 138 grad_berry(1,ipw,jband) = & 139 & grad_berry(1,ipw,jband) + & 140 & fac*qmat(1,iband,jband,ikpt,1,idir)*wfr - fac*qmat(2,iband,jband,ikpt,1,idir)*wfi 141 142 grad_berry(2,ipw,jband) = & 143 & grad_berry(2,ipw,jband) + & 144 & fac*qmat(1,iband,jband,ikpt,1,idir)*wfi + fac*qmat(2,iband,jband,ikpt,1,idir)*wfr 145 end do 146 end do 147 end if 148 end do 149 150 ! compute <u^(0)_{k_j}|u^(1)_{k_j+1}> matrix---------------------------------------------------- 151 152 ! prepare to calculate overlap matrix 153 ikptn = ikpt 154 ikpt1 = dtefield%ikpt_dk(ikpt,1,idir) 155 icg = dtefield%cgindex(ikptn,isppol) 156 icg1 = dtefield%cgindex(ikpt1,isppol+nsppol) 157 npw_k1 = npwarr(ikptn) 158 npw_k2 = npwar1(ikpt1) 159 pwind_tmp(1:npw_k1) = pwindall((ikptn-1)*mpw_tmp+1:(ikptn-1)*mpw_tmp+npw_k1,7,idir) 160 161 vect1(:,0) = zero ; vect2(:,0) = zero 162 do jband = 1, dtefield%mband_occ 163 vect2(:,1:npw_k2) = & 164 & cg1(:,icg1 + 1 + (jband-1)*npw_k2*nspinor:icg1 + jband*npw_k2*nspinor) 165 if (npw_k2 < mpw_tmp) vect2(:,npw_k2+1:mpw_tmp) = zero 166 do iband = 1, dtefield%mband_occ 167 pwmin = (iband-1)*npw_k1*nspinor 168 pwmax = pwmin + npw_k1*nspinor 169 vect1(:,1:npw_k1) = & 170 & cg(:,icg + 1 + pwmin:icg + pwmax) 171 if (npw_k1 < mpw_tmp) vect1(:,npw_k1+1:mpw_tmp) = zero 172 call overlap_g(doti,dotr,mpw_tmp,npw_k1,npw_k2,nspinor,pwind_tmp,& 173 & vect1,vect2) 174 s1mat(1,iband,jband) = dotr 175 s1mat(2,iband,jband) = doti 176 end do ! iband 177 end do !jband 178 179 ! compute <u^(1)_{k_j}|u^(0)_{k_j+1}> matrix----------------------------------------------------- 180 181 ! prepare to calculate overlap matrix 182 ikpt1 = dtefield%ikpt_dk(ikpt,1,idir) 183 icg = dtefield%cgindex(ikpt,isppol+nsppol) 184 icg1 = dtefield%cgindex(ikpt1,isppol) 185 npw_k1 = npwar1(ikpt) 186 npw_k2 = npwarr(ikpt1) 187 pwind_tmp(1:npw_k1) = pwindall((ikpt-1)*mpw_tmp+1:(ikpt-1)*mpw_tmp+npw_k1,5,idir) 188 189 vect1(:,0) = zero ; vect2(:,0) = zero 190 do jband = 1, dtefield%mband_occ 191 vect2(:,1:npw_k2) = & 192 & cg(:,icg1 + 1 + (jband-1)*npw_k2*nspinor:icg1 + jband*npw_k2*nspinor) 193 if (npw_k2 < mpw_tmp) vect2(:,npw_k2+1:mpw_tmp) = zero 194 do iband = 1, dtefield%mband_occ 195 pwmin = (iband-1)*npw_k1*nspinor 196 pwmax = pwmin + npw_k1*nspinor 197 vect1(:,1:npw_k1) = & 198 & cg1(:,icg + 1 + pwmin:icg + pwmax) 199 if (npw_k1 < mpw_tmp) vect1(:,npw_k1+1:mpw_tmp) = zero 200 call overlap_g(doti,dotr,mpw_tmp,npw_k1,npw_k2,nspinor,pwind_tmp,& 201 & vect1,vect2) 202 s1mat(1,jband,iband) = s1mat(1,jband,iband) + dotr 203 s1mat(2,jband,iband) = s1mat(2,jband,iband) + doti 204 end do ! iband 205 end do !jband 206 207 Amat(:,:,:)=zero 208 209 ! calculate Amat 210 do iband=1, dtefield%mband_occ 211 do jband=1, dtefield%mband_occ 212 do kband=1, dtefield%mband_occ 213 Amat(1,iband,jband) = Amat(1,iband,jband) + s1mat(1,iband,kband)*qmat(1,kband,jband,ikpt,1,idir)& 214 & - s1mat(2,iband,kband)*qmat(2,kband,jband,ikpt,1,idir) 215 Amat(2,iband,jband) = Amat(2,iband,jband) + s1mat(1,iband,kband)*qmat(2,kband,jband,ikpt,1,idir)& 216 & + s1mat(2,iband,kband)*qmat(1,kband,jband,ikpt,1,idir) 217 end do 218 end do 219 end do 220 221 Bmat(:,:,:)=zero 222 223 ! calculate Bmat 224 ikptn = dtefield%ikpt_dk(ikpt,7,idir) 225 do iband=1, dtefield%mband_occ 226 do jband=1, dtefield%mband_occ 227 do kband=1, dtefield%mband_occ 228 Bmat(1,jband,kband) = Bmat(1,jband,kband) + Amat(1,iband,kband)*qmat(1,jband,iband,ikptn,1,idir)& 229 & - Amat(2,iband,kband)*qmat(2,jband,iband,ikptn,1,idir) 230 Bmat(2,jband,kband) = Bmat(2,jband,kband) + Amat(1,iband,kband)*qmat(2,jband,iband,ikptn,1,idir)& 231 & + Amat(2,iband,kband)*qmat(1,jband,iband,ikptn,1,idir) 232 end do 233 end do 234 end do 235 236 ! calc. the second term of gradient------------------------------ 237 238 ! preparation 239 240 ikptnp1 = dtefield%ikpt_dk(ikpt,3,idir) 241 icg = dtefield%cgindex(ikptnp1,isppol) 242 npw_k1 = npwar1(ikpt) 243 npw_k2 = npwarr(ikptnp1) 244 pwind_tmp(1:npw_k1) = pwindall((ikpt-1)*mpw_tmp+1:(ikpt-1)*mpw_tmp+npw_k1,5,idir) 245 246 z1(:) = zero 247 z2(:) = zero 248 249 do ipw = 1, npw_k1 250 jpw = pwind_tmp(ipw) 251 if (jpw > 0) then 252 do iband = 1, dtefield%mband_occ 253 wfr = cg(1,icg + (iband - 1)*npw_k2*nspinor + jpw) 254 wfi = cg(2,icg + (iband - 1)*npw_k2*nspinor + jpw) 255 do jband=1, dtefield%mband_occ 256 grad_berry(1,ipw,jband) = grad_berry(1,ipw,jband) - fac*(Bmat(1,iband,jband)*wfr - Bmat(2,iband,jband)*wfi) 257 grad_berry(2,ipw,jband) = grad_berry(2,ipw,jband) - fac*(Bmat(1,iband,jband)*wfi + Bmat(2,iband,jband)*wfr) 258 end do 259 end do 260 end if 261 end do 262 263 ! Second part of gradient of Berry phase++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 264 265 vect1(:,0) = zero ; vect2(:,0) = zero 266 267 ! prepare 268 ikpt1 = dtefield%ikpt_dk(ikpt,2,idir) 269 icg1 = dtefield%cgindex(ikpt1,isppol+nsppol) 270 npw_k1 = npwar1(ikpt) 271 npw_k2 = npwar1(ikpt1) 272 pwind_tmp(1:npw_k1) = pwindall((ikpt-1)*mpw_tmp+1:(ikpt-1)*mpw_tmp+npw_k1,4,idir) 273 274 do ipw = 1, npw_k1 275 jpw = pwind_tmp(ipw) 276 if (jpw > 0) then 277 do iband = 1, dtefield%mband_occ 278 wfr = cg1(1,icg1 + (iband - 1)*npw_k2*nspinor + jpw) 279 wfi = cg1(2,icg1 + (iband - 1)*npw_k2*nspinor + jpw) 280 do jband = 1, dtefield%mband_occ 281 grad_berry(1,ipw,jband) = & 282 & grad_berry(1,ipw,jband) - & 283 & fac*qmat(1,iband,jband,ikpt,2,idir)*wfr + fac*qmat(2,iband,jband,ikpt,2,idir)*wfi 284 285 grad_berry(2,ipw,jband) = & 286 & grad_berry(2,ipw,jband) - & 287 & fac*qmat(1,iband,jband,ikpt,2,idir)*wfi - fac*qmat(2,iband,jband,ikpt,2,idir)*wfr 288 end do 289 end do 290 end if 291 end do 292 293 ! compute <u^(0)_{k_j}|u^(1)_{k_j-1}> matrix---------------------------------------------------- 294 295 ! prepare to calculate overlap matrix 296 ikptn = ikpt 297 ikpt1 = dtefield%ikpt_dk(ikpt,2,idir) 298 icg = dtefield%cgindex(ikptn,isppol) 299 icg1 = dtefield%cgindex(ikpt1,isppol+nsppol) 300 npw_k1 = npwarr(ikptn) 301 npw_k2 = npwar1(ikpt1) 302 pwind_tmp(1:npw_k1) =pwindall((ikptn-1)*mpw_tmp+1:(ikptn-1)*mpw_tmp+npw_k1,8,idir) 303 304 vect1(:,0) = zero ; vect2(:,0) = zero 305 do jband = 1, dtefield%mband_occ 306 vect2(:,1:npw_k2) = & 307 & cg1(:,icg1 + 1 + (jband-1)*npw_k2*nspinor:icg1 + jband*npw_k2*nspinor) 308 if (npw_k2 < mpw_tmp) vect2(:,npw_k2+1:mpw_tmp) = zero 309 do iband = 1, dtefield%mband_occ 310 pwmin = (iband-1)*npw_k1*nspinor 311 pwmax = pwmin + npw_k1*nspinor 312 vect1(:,1:npw_k1) = & 313 & cg(:,icg + 1 + pwmin:icg + pwmax) 314 if (npw_k1 < mpw_tmp) vect1(:,npw_k1+1:mpw_tmp) = zero 315 call overlap_g(doti,dotr,mpw_tmp,npw_k1,npw_k2,nspinor,pwind_tmp,& 316 & vect1,vect2) 317 s1mat(1,iband,jband) = dotr 318 s1mat(2,iband,jband) = doti 319 end do ! iband 320 end do !jband 321 322 ! compute <u^(1)_{k_j}|u^(0)_{k_j-1}> matrix----------------------------------------------------- 323 324 ! prepare to calculate overlap matrix 325 ikpt1 = dtefield%ikpt_dk(ikpt,2,idir) 326 icg = dtefield%cgindex(ikpt,isppol+nsppol) 327 icg1 = dtefield%cgindex(ikpt1,isppol) 328 npw_k1 = npwarr(ikpt) 329 npw_k2 = npwar1(ikpt1) 330 pwind_tmp(1:npw_k1) =pwindall((ikpt-1)*mpw_tmp+1:(ikpt-1)*mpw_tmp+npw_k1,6,idir) 331 332 vect1(:,0) = zero ; vect2(:,0) = zero 333 do jband = 1, dtefield%mband_occ 334 vect2(:,1:npw_k2) = & 335 & cg(:,icg1 + 1 + (jband-1)*npw_k2*nspinor:icg1 + jband*npw_k2*nspinor) 336 if (npw_k2 < mpw_tmp) vect2(:,npw_k2+1:mpw_tmp) = zero 337 do iband = 1, dtefield%mband_occ 338 pwmin = (iband-1)*npw_k1*nspinor 339 pwmax = pwmin + npw_k1*nspinor 340 vect1(:,1:npw_k1) = & 341 & cg1(:,icg + 1 + pwmin:icg + pwmax) 342 if (npw_k1 < mpw_tmp) vect1(:,npw_k1+1:mpw_tmp) = zero 343 call overlap_g(doti,dotr,mpw_tmp,npw_k1,npw_k2,nspinor,pwind_tmp,& 344 & vect1,vect2) 345 s1mat(1,jband,iband) = s1mat(1,jband,iband) + dotr 346 s1mat(2,jband,iband) = s1mat(2,jband,iband) + doti 347 end do ! iband 348 end do !jband 349 350 Amat(:,:,:)=zero 351 352 ! calculate Amat 353 do iband=1, dtefield%mband_occ 354 do jband=1, dtefield%mband_occ 355 do kband=1, dtefield%mband_occ 356 Amat(1,iband,jband) = Amat(1,iband,jband) + s1mat(1,iband,kband)*qmat(1,kband,jband,ikpt,2,idir)& 357 & - s1mat(2,iband,kband)*qmat(2,kband,jband,ikpt,2,idir) 358 Amat(2,iband,jband) = Amat(2,iband,jband) + s1mat(1,iband,kband)*qmat(2,kband,jband,ikpt,2,idir)& 359 & + s1mat(2,iband,kband)*qmat(1,kband,jband,ikpt,2,idir) 360 end do 361 end do 362 end do 363 364 Bmat(:,:,:)=zero 365 366 ! calculate Bmat 367 ikptn = dtefield%ikpt_dk(ikpt,7,idir) 368 do iband=1, dtefield%mband_occ 369 do jband=1, dtefield%mband_occ 370 do kband=1, dtefield%mband_occ 371 Bmat(1,jband,kband) = Bmat(1,jband,kband) + Amat(1,iband,kband)*qmat(1,jband,iband,ikptn,2,idir)& 372 & - Amat(2,iband,kband)*qmat(2,jband,iband,ikptn,2,idir) 373 Bmat(2,jband,kband) = Bmat(2,jband,kband) + Amat(1,iband,kband)*qmat(2,jband,iband,ikptn,2,idir)& 374 + Amat(2,iband,kband)*qmat(1,jband,iband,ikptn,2,idir) 375 end do 376 end do 377 end do 378 379 ! calc. the second term of gradient------------------------------ 380 381 ! preparation 382 383 ikptnp1 = dtefield%ikpt_dk(ikpt,4,idir) 384 icg = dtefield%cgindex(ikptnp1,isppol) 385 npw_k1 = npwar1(ikpt) 386 npw_k2 = npwarr(ikptnp1) 387 pwind_tmp(1:npw_k1) =pwindall((ikpt-1)*mpw_tmp+1:(ikpt-1)*mpw_tmp+npw_k1,6,idir) 388 z1(:) = zero 389 z2(:) = zero 390 do ipw = 1, npw_k1 391 jpw = pwind_tmp(ipw) 392 if (jpw > 0) then 393 do iband = 1, dtefield%mband_occ 394 wfr = cg(1,icg + (iband - 1)*npw_k2*nspinor + jpw) 395 wfi = cg(2,icg + (iband - 1)*npw_k2*nspinor + jpw) 396 do jband=1, dtefield%mband_occ 397 grad_berry(1,ipw,jband) = grad_berry(1,ipw,jband) + fac*(Bmat(1,iband,jband)*wfr - Bmat(2,iband,jband)*wfi) 398 grad_berry(2,ipw,jband) = grad_berry(2,ipw,jband) + fac*(Bmat(1,iband,jband)*wfi + Bmat(2,iband,jband)*wfr) 399 end do 400 end do 401 end if 402 end do 403 404 end do !idir 405 406 !!----------------------------------------third part of gradient------------------------------------------------------ 407 do idir=1,3 408 fac = rprimd(idir_efield,idir)*dble(nkpt)/& 409 & (dble(dtefield%nstr(idir))*four_pi) 410 411 ! prepare 412 ikpt1 = dtefield%ikpt_dk(ikpt,1,idir) 413 icg1 = dtefield%cgindex(ikpt1,isppol) 414 npw_k1 = npwar1(ikpt) 415 npw_k2 = npwarr(ikpt1) 416 pwind_tmp(1:npw_k1) = pwindall((ikpt-1)*mpw_tmp+1:(ikpt-1)*mpw_tmp+npw_k1,5,idir) 417 do ipw = 1, npw_k1 418 jpw = pwind_tmp(ipw) 419 if (jpw > 0) then 420 do iband = 1, dtefield%mband_occ 421 wfr = cg(1,icg1 + (iband - 1)*npw_k2*nspinor + jpw) 422 wfi = cg(2,icg1 + (iband - 1)*npw_k2*nspinor + jpw) 423 do jband = 1, dtefield%mband_occ 424 grad_berry(1,ipw,jband) = & 425 & grad_berry(1,ipw,jband) + & 426 & fac*qmat(1,iband,jband,ikpt,1,idir)*wfr - fac*qmat(2,iband,jband,ikpt,1,idir)*wfi 427 grad_berry(2,ipw,jband) = & 428 & grad_berry(2,ipw,jband) + & 429 & fac*qmat(1,iband,jband,ikpt,1,idir)*wfi + fac*qmat(2,iband,jband,ikpt,1,idir)*wfr 430 end do 431 end do 432 end if 433 end do 434 435 ! prepare 436 ikpt1 = dtefield%ikpt_dk(ikpt,2,idir) 437 icg1 = dtefield%cgindex(ikpt1,isppol) 438 npw_k1 = npwar1(ikpt) 439 npw_k2 = npwarr(ikpt1) 440 pwind_tmp(1:npw_k1) = pwindall((ikpt-1)*mpw_tmp+1:(ikpt-1)*mpw_tmp+npw_k1,6,idir) 441 442 do ipw = 1, npw_k1 443 jpw = pwind_tmp(ipw) 444 if (jpw > 0) then 445 do iband = 1, dtefield%mband_occ 446 wfr = cg(1,icg1 + (iband - 1)*npw_k2*nspinor + jpw) 447 wfi = cg(2,icg1 + (iband - 1)*npw_k2*nspinor + jpw) 448 do jband = 1, dtefield%mband_occ 449 grad_berry(1,ipw,jband) = & 450 & grad_berry(1,ipw,jband) - & 451 & fac*qmat(1,iband,jband,ikpt,2,idir)*wfr + fac*qmat(2,iband,jband,ikpt,2,idir)*wfi 452 453 grad_berry(2,ipw,jband) = & 454 & grad_berry(2,ipw,jband) - & 455 & fac*qmat(1,iband,jband,ikpt,2,idir)*wfi - fac*qmat(2,iband,jband,ikpt,2,idir)*wfr 456 457 end do 458 end do 459 end if 460 end do 461 462 end do !idir 463 464 ABI_DEALLOCATE(vect1) 465 ABI_DEALLOCATE(vect2) 466 ABI_DEALLOCATE(s1mat) 467 ABI_DEALLOCATE(Amat) 468 ABI_DEALLOCATE(Bmat) 469 ABI_DEALLOCATE(pwind_tmp) 470 471 end subroutine dfptff_gbefd