TABLE OF CONTENTS
ABINIT/dfptff_die [ Functions ]
NAME
dfptff_die
FUNCTION
calculate electric susceptibility tensor in Eq.(28) 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 idirpert = the current coloumn of the dielectric permittivity tensor 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> rprimd(3,3)=dimensional primitive translations in real space (bohr)
OUTPUT
diet = electric susceptibility tensor
PARENTS
dfpt_scfcv
CHILDREN
overlap_g
SOURCE
55 #if defined HAVE_CONFIG_H 56 #include "config.h" 57 #endif 58 59 #include "abi_common.h" 60 61 subroutine dfptff_die(cg,cg1,dtefield,d2lo,idirpert,ipert,mband,mkmem,& 62 & mpw,mpw1,mpert,nkpt,npwarr,npwar1,nsppol,nspinor,pwindall,qmat,rprimd) 63 64 use defs_basis 65 use m_profiling_abi 66 use m_efield 67 68 use m_cgtools, only : overlap_g 69 70 !This section has been created automatically by the script Abilint (TD). 71 !Do not modify the following lines by hand. 72 #undef ABI_FUNC 73 #define ABI_FUNC 'dfptff_die' 74 !End of the abilint section 75 76 implicit none 77 78 !Arguments ---------------------------------------- 79 !scalars 80 integer,intent(in) :: idirpert,ipert,mband,mkmem,mpert,mpw,mpw1,nkpt,nspinor 81 integer,intent(in) :: nsppol 82 type(efield_type),intent(in) :: dtefield 83 !arrays 84 integer,intent(in) :: npwar1(nkpt),npwarr(nkpt) 85 integer,intent(in) :: pwindall(max(mpw,mpw1)*mkmem,8,3) 86 real(dp),intent(in) :: cg(2,mpw*mband*mkmem*nspinor*nsppol) 87 real(dp),intent(in) :: cg1(2,mpw1*mband*mkmem*nspinor*nsppol) 88 real(dp),intent(in) :: qmat(2,dtefield%mband_occ,dtefield%mband_occ,nkpt,2,3) 89 real(dp),intent(in) :: rprimd(3,3) 90 real(dp),intent(inout) :: d2lo(2,3,mpert,3,mpert) !vz_i 91 92 !Local variables ---------------------------------- 93 !scalars 94 integer :: ialpha,iband,icg,icg1,idir,ikpt,ikpt1,jband,mpw_tmp,npw_k1 95 integer :: npw_k2,pwmax,pwmin 96 real(dp) :: doti,dotr,e0,fac 97 !arrays 98 integer,allocatable :: pwind_tmp(:) 99 real(dp) :: edir(3) 100 real(dp),allocatable :: s1mat(:,:,:),vect1(:,:),vect2(:,:) 101 102 ! ************************************************************************* 103 104 !calculate s1 matrices ----------------------------- 105 mpw_tmp=max(mpw,mpw1) 106 ABI_ALLOCATE(s1mat,(2,dtefield%mband_occ,dtefield%mband_occ)) 107 ABI_ALLOCATE(vect1,(2,0:mpw_tmp)) 108 ABI_ALLOCATE(vect2,(2,0:mpw_tmp)) 109 ABI_ALLOCATE(pwind_tmp,(mpw_tmp)) 110 vect1(:,0) = zero ; vect2(:,0) = zero 111 112 edir(:)=zero 113 114 do ikpt=1,nkpt 115 do idir=1,3 116 ! compute <u^(0)_{k_j}|u^(1)_{k_j+1,q}> matrix--- q=0 ---------------------------------------- 117 118 ! prepare to calculate overlap matrix 119 ikpt1 = dtefield%ikpt_dk(ikpt,1,idir) 120 icg = dtefield%cgindex(ikpt,1) 121 icg1 = dtefield%cgindex(ikpt1,1+nsppol) 122 npw_k1 = npwarr(ikpt) 123 npw_k2 = npwar1(ikpt1) 124 pwind_tmp(1:npw_k1) = pwindall((ikpt-1)*mpw_tmp+1:(ikpt-1)*mpw_tmp+npw_k1,7,idir) 125 126 vect1(:,0) = zero ; vect2(:,0) = zero 127 do jband = 1, dtefield%mband_occ 128 vect2(:,1:npw_k2) = & 129 & cg1(:,icg1 + 1 + (jband-1)*npw_k2*nspinor:icg1 + jband*npw_k2*nspinor) 130 if (npw_k2 < mpw_tmp) vect2(:,npw_k2+1:mpw_tmp) = zero 131 do iband = 1, dtefield%mband_occ 132 pwmin = (iband-1)*npw_k1*nspinor 133 pwmax = pwmin + npw_k1*nspinor 134 vect1(:,1:npw_k1) = & 135 & cg(:,icg + 1 + pwmin:icg + pwmax) 136 if (npw_k1 < mpw_tmp) vect1(:,npw_k1+1:mpw_tmp) = zero 137 call overlap_g(doti,dotr,mpw_tmp,npw_k1,npw_k2,nspinor,pwind_tmp,& 138 & vect1,vect2) 139 s1mat(1,iband,jband) = dotr 140 s1mat(2,iband,jband) = doti 141 end do ! iband 142 end do !jband 143 144 ! compute <u^(1)_{k_j,q}|u^(0)_{k_j+1}> matrix-- q=0 ------------------------------------- 145 146 ! prepare to calculate overlap matrix 147 ikpt1 = dtefield%ikpt_dk(ikpt,1,idir) 148 icg = dtefield%cgindex(ikpt,1+nsppol) 149 icg1 = dtefield%cgindex(ikpt1,1) 150 npw_k1 = npwar1(ikpt) 151 npw_k2 = npwarr(ikpt1) 152 pwind_tmp(1:npw_k1) = pwindall((ikpt-1)*mpw_tmp+1:(ikpt-1)*mpw_tmp+npw_k1,5,idir) 153 vect1(:,0) = zero ; vect2(:,0) = zero 154 do jband = 1, dtefield%mband_occ 155 vect2(:,1:npw_k2) = & 156 & cg(:,icg1 + 1 + (jband-1)*npw_k2*nspinor:icg1 + jband*npw_k2*nspinor) 157 if (npw_k2 < mpw_tmp) vect2(:,npw_k2+1:mpw_tmp) = zero 158 do iband = 1, dtefield%mband_occ 159 pwmin = (iband-1)*npw_k1*nspinor 160 pwmax = pwmin + npw_k1*nspinor 161 vect1(:,1:npw_k1) = & 162 & cg1(:,icg + 1 + pwmin:icg + pwmax) 163 if (npw_k1 < mpw_tmp) vect1(:,npw_k1+1:mpw_tmp) = zero 164 call overlap_g(doti,dotr,mpw_tmp,npw_k1,npw_k2,nspinor,pwind_tmp,& 165 & vect1,vect2) 166 s1mat(1,iband,jband) = s1mat(1,iband,jband) + dotr 167 s1mat(2,iband,jband) = s1mat(2,iband,jband) + doti 168 end do ! iband 169 end do !jband 170 171 ! sum over the whole------------------------------------------------------------ 172 173 e0=zero 174 175 do iband=1,dtefield%mband_occ 176 do jband=1,dtefield%mband_occ 177 e0 = e0 + (s1mat(1,iband,jband)*qmat(2,jband,iband,ikpt,1,idir)& 178 & + s1mat(2,iband,jband)*qmat(1,jband,iband,ikpt,1,idir)) 179 180 end do 181 end do 182 183 do ialpha=1,3 184 fac = rprimd(ialpha,idir)/& 185 & (dble(dtefield%nstr(idir))*pi) 186 edir(ialpha)=edir(ialpha)+ e0*fac 187 end do 188 189 end do !idir 190 end do !ikpt 191 192 d2lo(1,1:3,ipert,idirpert,ipert)=edir(:) 193 194 ABI_DEALLOCATE(s1mat) 195 ABI_DEALLOCATE(vect1) 196 ABI_DEALLOCATE(vect2) 197 ABI_DEALLOCATE(pwind_tmp) 198 end subroutine dfptff_die