TABLE OF CONTENTS
ABINIT/calc_optical_mels [ Functions ]
NAME
calc_optical_mels
FUNCTION
Calculate all optical matrix elements in the BZ.
INPUTS
lomo_spin(Wfd%nsppol)=Index of the lomo band for the different spins. lomo_min,max_band=minimum and max band index to be calculated. nkbz=Number of points in the full Brillouin zone. inclvkb=if different from 0, [Vnl,r] is included in the calculation of the matrix element of the velocity operator. No meaning for PAW (except for DFT+U) qpt(3) Kmesh<kmesh_t>=Info on the k-point sampling for wave functions. Cryst<crystal_t>=Structure defining the crystalline structure. KS_Bst<ebands_t> Pawtab(Cryst%ntypat*usepaw)<pawtab_type>=PAW tabulated starting data Psps <pseudopotential_type>=variables related to pseudopotentials. Hur(Cryst%natom*usepaw)<pawhur_t>=Only for PAW and DFT+U, quantities used to evaluate the commutator [H_u,r]. Wfd<wfdgw_t>=Handler for the wavefunctions.
OUTPUT
opt_cvk(lomo_min:max_band,lomo_min:max_band,nkbz,nsppol)=Matrix elements <c k|e^{+iqr}|v k>
SOURCE
77 subroutine calc_optical_mels(Wfd,Kmesh,KS_Bst,Cryst,Psps,Pawtab,Hur,& 78 & inclvkb,lomo_spin,lomo_min,max_band,nkbz,qpoint,opt_cvk) 79 80 !Arguments ------------------------------------ 81 !scalars 82 integer,intent(in) :: nkbz,inclvkb,lomo_min,max_band 83 type(kmesh_t),intent(in) :: Kmesh 84 type(crystal_t),intent(in) :: Cryst 85 type(pseudopotential_type),intent(in) :: Psps 86 type(ebands_t),intent(in) :: KS_Bst 87 type(wfdgw_t),target,intent(inout) :: Wfd 88 !arrays 89 integer,intent(in) :: lomo_spin(Wfd%nsppol) 90 real(dp),intent(in) :: qpoint(3) 91 complex(dpc),intent(out) :: opt_cvk(lomo_min:max_band,lomo_min:max_band,nkbz,Wfd%nsppol) 92 type(pawtab_type),intent(in) :: Pawtab(Cryst%ntypat*Wfd%usepaw) 93 type(pawhur_t),intent(in) :: Hur(Cryst%natom*Wfd%usepaw) 94 95 !Local variables ------------------------------ 96 !scalars 97 integer :: nsppol,usepaw,nspinor,comm,spin,npw_k,istwf_k,my_nbbp 98 integer :: ik_bz,ik_ibz,itim_k,isym_k,ib_c,ib_v,ierr,my_rank 99 real(dp) :: ediff 100 complex(dpc) :: emcvk 101 character(len=500) :: msg 102 type(vkbr_t) :: vkbr 103 type(wave_t),pointer :: wave_v, wave_c 104 !arrays 105 integer,allocatable :: bbp_distrb(:,:) 106 integer,ABI_CONTIGUOUS pointer :: kg_k(:,:) 107 real(dp) :: mat_dp(3,3),qrot(3),b1(3),b2(3),b3(3),kbz(3) 108 complex(dpc),allocatable :: ir_kibz(:,:,:,:,:) 109 complex(gwpc), ABI_CONTIGUOUS pointer :: ug_c(:),ug_v(:) 110 complex(gwpc) :: ihrc(3,Wfd%nspinor**2) 111 logical :: bbp_mask(Wfd%mband,Wfd%mband) 112 type(pawcprj_type),allocatable :: Cp_v(:,:),Cp_c(:,:) 113 114 !************************************************************************ 115 116 call wrtout(std_out," Calculating optical matrix elements in the IBZ","COLL") 117 ABI_CHECK(Wfd%nspinor==1,"nspinor==2 not coded") 118 119 comm = Wfd%comm 120 my_rank = Wfd%my_rank 121 122 nsppol = Wfd%nsppol 123 nspinor = Wfd%nspinor 124 usepaw = Wfd%usepaw 125 126 if (usepaw==1) then 127 ABI_MALLOC(Cp_v,(Wfd%natom,nspinor)) 128 call pawcprj_alloc(Cp_v,0,Wfd%nlmn_atm) 129 ABI_MALLOC(Cp_c,(Wfd%natom,nspinor)) 130 call pawcprj_alloc(Cp_c,0,Wfd%nlmn_atm) 131 end if 132 133 if (inclvkb==1.and.usepaw==0) then 134 ABI_ERROR("inclvkb==1 not coded,using inclvkb==2") 135 end if 136 ! 137 ! Calculate the matrix elements of ir in the IBZ. 138 ABI_MALLOC(ir_kibz,(3,lomo_min:max_band,lomo_min:max_band,Wfd%nkibz,nsppol)) 139 ir_kibz=czero 140 141 ABI_MALLOC(bbp_distrb, (Wfd%mband,Wfd%mband)) 142 143 do spin=1,nsppol 144 do ik_ibz=1,Wfd%nkibz 145 ! 146 ! Distribute the (b,b') entries. 147 bbp_mask=.FALSE.; bbp_mask(lomo_spin(spin):max_band,lomo_spin(spin):max_band)=.TRUE. 148 call wfd%distribute_bbp(ik_ibz,spin,"All",my_nbbp,bbp_distrb,bbp_mask=bbp_mask) 149 if (ALL(bbp_distrb/=my_rank)) CYCLE 150 151 istwf_k = Wfd%istwfk(ik_ibz) 152 ABI_CHECK(istwf_k==1,"istwf_k/=1 not coded") ! KB stuff is missing. 153 npw_k = Wfd%npwarr(ik_ibz) 154 kg_k => Wfd%Kdata(ik_ibz)%kg_k 155 156 if (inclvkb/=0.and.usepaw==0) then 157 ! Prepare term i <n,k|[Vnl,r]|n"k> 158 call vkbr_init(vkbr,Cryst,Psps,inclvkb,istwf_k,npw_k,Kmesh%ibz(:,ik_ibz),kg_k) 159 end if 160 161 ! Note: spinorial case is not coded therefore we work with ihrc(:,1). 162 ! TODO: The lower triangle can be Reconstructed by symmetry. 163 do ib_v=lomo_spin(spin),max_band ! Loop over bands 164 if ( ALL(bbp_distrb(ib_v,:)/=my_rank) ) CYCLE 165 166 ABI_CHECK(wfd%get_wave_ptr(ib_v, ik_ibz, spin, wave_v, msg) == 0, msg) 167 ug_v => wave_v%ug 168 if (usepaw==1) call wfd%get_cprj(ib_v,ik_ibz,spin,Cryst,Cp_v,sorted=.FALSE.) 169 170 do ib_c=lomo_spin(spin),max_band 171 if (bbp_distrb(ib_v,ib_c)/=my_rank) CYCLE 172 ABI_CHECK(wfd%get_wave_ptr(ib_c, ik_ibz, spin, wave_c, msg) == 0, msg) 173 ug_c => wave_c%ug 174 175 if (usepaw==0) then 176 ! Calculate matrix elements of i[H,r] for NC pseudopotentials. 177 ihrc = nc_ihr_comm(vkbr,cryst,psps,npw_k,nspinor,istwf_k,inclvkb,Kmesh%ibz(:,ik_ibz),ug_c,ug_v,kg_k) 178 179 else 180 ! Matrix elements of i[H,r] for PAW. 181 call wfd%get_cprj(ib_c,ik_ibz,spin,Cryst,Cp_c,sorted=.FALSE.) 182 183 ihrc = paw_ihr(spin,nspinor,npw_k,istwf_k,Kmesh%ibz(:,ik_ibz),Cryst,Pawtab,ug_c,ug_v,kg_k,Cp_c,Cp_v,HUr) 184 end if 185 ! 186 ! Save matrix elements of i*r in the IBZ 187 ediff = KS_Bst%eig(ib_c,ik_ibz,spin) - KS_BSt%eig(ib_v,ik_ibz,spin) 188 if (ABS(ediff)<tol16) ediff=tol6 ! Treat a possible degeneracy between v and c. 189 ir_kibz(:,ib_c,ib_v,ik_ibz,spin) = ihrc(:,1)/ediff 190 191 end do !ib_c 192 end do !ib_v 193 194 call vkbr_free(vkbr) 195 end do !spin 196 end do !ik_ibz 197 198 ! Collect results on each node. 199 call xmpi_sum(ir_kibz,comm,ierr) 200 201 ABI_FREE(bbp_distrb) 202 203 if (usepaw==1) then 204 call pawcprj_free(Cp_v) 205 ABI_FREE(Cp_v) 206 call pawcprj_free(Cp_c) 207 ABI_FREE(Cp_c) 208 end if 209 ! 210 ! ====================================================== 211 ! ==== Calculate Fcv(kBZ) in the full Brilouin zone ==== 212 ! ====================================================== 213 ! 214 ! Symmetrization of the matrix elements of the position operator. 215 ! <Sk b|r|Sk b'> = R <k b|r|k b'> + \tau \delta_{bb'} 216 ! where S is one of the symrec operations in reciprocal space, R is the 217 ! corresponding operation in real space, \tau being the associated fractional translations. 218 ! 219 ! q.Mcv( Sk) = S^{-1}q. Mcv(k) 220 ! q.Mcv(-Sk) = -S^{-1}q. CONJG(Mcv(k)) if time-reversal is used. 221 222 b1=Cryst%gprimd(:,1)*two_pi 223 b2=Cryst%gprimd(:,2)*two_pi 224 b3=Cryst%gprimd(:,3)*two_pi 225 226 opt_cvk = czero 227 do spin=1,nsppol 228 do ik_bz=1,nkbz 229 ! 230 ! Get ik_ibz, and symmetries index from ik_bz. 231 call Kmesh%get_BZ_item(ik_bz,kbz,ik_ibz,isym_k,itim_k) 232 233 mat_dp = DBLE(Cryst%symrec(:,:,isym_k)) 234 call matrginv(mat_dp,3,3) ! Invert 235 qrot = (3-2*itim_k) * MATMUL(mat_dp,qpoint) 236 237 do ib_v=lomo_spin(spin),max_band ! Loops over the bands C and V start 238 do ib_c=lomo_spin(spin),max_band 239 !if (ib_c==ib_v) CYCLE 240 emcvk = pdtqrc(qrot,ir_kibz(:,ib_c,ib_v,ik_ibz,spin),b1,b2,b3) 241 if (itim_k==2) emcvk = CONJG(emcvk) 242 opt_cvk(ib_c,ib_v,ik_bz,spin) = emcvk 243 end do !ib_c 244 end do !ib_v 245 246 end do !ik_bz 247 end do !spin 248 249 ABI_FREE(ir_kibz) 250 251 call xmpi_barrier(comm) 252 253 contains
ABINIT/m_wfd_optic [ Modules ]
NAME
m_wfd_optic
FUNCTION
Functions to compute optical matrix elements using the wavefunction descriptor.
COPYRIGHT
Copyright (C) 2008-2022 ABINIT group (MG) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
SOURCE
16 #if defined HAVE_CONFIG_H 17 #include "config.h" 18 #endif 19 20 #include "abi_common.h" 21 22 module m_wfd_optic 23 24 use defs_basis 25 use m_errors 26 use m_abicore 27 use m_xmpi 28 29 use defs_datatypes, only : ebands_t, pseudopotential_type 30 use m_hide_lapack, only : matrginv 31 use m_bz_mesh, only : kmesh_t 32 use m_crystal, only : crystal_t 33 use m_vkbr, only : vkbr_t, vkbr_free, vkbr_init, nc_ihr_comm 34 use m_wfd, only : wfdgw_t, wave_t 35 use m_pawtab, only : pawtab_type 36 use m_pawcprj, only : pawcprj_type, pawcprj_alloc, pawcprj_free 37 use m_paw_hr, only : pawhur_t, paw_ihr 38 39 implicit none 40 41 private
ABINIT/pdtqrc [ Functions ]
NAME
pdtqrc
FUNCTION
Calculate the dot product of a real vector with a complex vector, where each is in terms of b1-b3
INPUTS
OUTPUT
SOURCE
269 pure function pdtqrc(R,C,b1,b2,b3) 270 271 !Arguments ------------------------------------ 272 !arrays 273 real(dp),intent(in) :: R(3),b1(3),b2(3),b3(3) 274 complex(dpc),intent(in) :: C(3) 275 complex(dpc) :: pdtqrc 276 277 !Local variables ------------------------------ 278 !scalars 279 integer :: ii 280 281 !************************************************************************ 282 283 pdtqrc=czero 284 do ii=1,3 285 pdtqrc = pdtqrc + (R(1)*b1(ii)+R(2)*b2(ii)+R(3)*b3(ii)) * & 286 & (C(1)*b1(ii)+C(2)*b2(ii)+C(3)*b3(ii)) 287 end do 288 289 end function pdtqrc