TABLE OF CONTENTS


ABINIT/calc_optical_mels [ Functions ]

[ Top ] [ 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 LDA+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 LDA+U, quantities used to evaluate the commutator [H_u,r].
 Wfd<wfd_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>

PARENTS

      m_exc_spectra,m_haydock

CHILDREN

      get_bz_item,matrginv,pawcprj_alloc,pawcprj_free,vkbr_free,vkbr_init
      wfd_distribute_bbp,wfd_get_cprj,wrtout,xmpi_barrier,xmpi_sum

SOURCE

 89 subroutine calc_optical_mels(Wfd,Kmesh,KS_Bst,Cryst,Psps,Pawtab,Hur,&
 90 &  inclvkb,lomo_spin,lomo_min,max_band,nkbz,qpoint,opt_cvk)
 91 
 92 
 93 !This section has been created automatically by the script Abilint (TD).
 94 !Do not modify the following lines by hand.
 95 #undef ABI_FUNC
 96 #define ABI_FUNC 'calc_optical_mels'
 97 !End of the abilint section
 98 
 99  implicit none
100 
101 !Arguments ------------------------------------
102 !scalars
103  integer,intent(in) :: nkbz,inclvkb,lomo_min,max_band
104  type(kmesh_t),intent(in) :: Kmesh
105  type(crystal_t),intent(in) :: Cryst
106  type(pseudopotential_type),intent(in) :: Psps
107  type(ebands_t),intent(in) :: KS_Bst
108  type(wfd_t),target,intent(inout) :: Wfd
109 !arrays
110  integer,intent(in) :: lomo_spin(Wfd%nsppol)
111  real(dp),intent(in) :: qpoint(3)
112  complex(dpc),intent(out) :: opt_cvk(lomo_min:max_band,lomo_min:max_band,nkbz,Wfd%nsppol)
113  type(pawtab_type),intent(in) :: Pawtab(Cryst%ntypat*Wfd%usepaw)
114  type(pawhur_t),intent(in) :: Hur(Cryst%natom*Wfd%usepaw)
115 
116 !Local variables ------------------------------
117 !scalars
118  integer :: nsppol,usepaw,nspinor,comm,spin,npw_k,istwf_k,my_nbbp
119  integer :: ik_bz,ik_ibz,itim_k,isym_k,ib_c,ib_v,ierr,my_rank
120  real(dp) :: ediff
121  complex(dpc) :: emcvk
122  type(vkbr_t) :: vkbr
123 !arrays
124  integer,allocatable :: bbp_distrb(:,:)
125  integer,ABI_CONTIGUOUS pointer :: kg_k(:,:)
126  real(dp) :: mat_dp(3,3),qrot(3),b1(3),b2(3),b3(3),kbz(3)
127  complex(dpc),allocatable :: ir_kibz(:,:,:,:,:)
128  complex(gwpc), ABI_CONTIGUOUS pointer :: ug_c(:),ug_v(:)
129  complex(gwpc) :: ihrc(3,Wfd%nspinor**2)
130  logical :: bbp_mask(Wfd%mband,Wfd%mband)
131  type(pawcprj_type),allocatable :: Cp_v(:,:),Cp_c(:,:)
132 
133 !************************************************************************
134 
135  call wrtout(std_out," Calculating optical matrix elements in the IBZ","COLL")
136  ABI_CHECK(Wfd%nspinor==1,"nspinor==2 not coded")
137 
138  comm = Wfd%comm
139  my_rank = Wfd%my_rank
140 
141  nsppol  = Wfd%nsppol
142  nspinor = Wfd%nspinor
143  usepaw  = Wfd%usepaw
144 
145  if (usepaw==1) then
146    ABI_DT_MALLOC(Cp_v,(Wfd%natom,nspinor))
147    call pawcprj_alloc(Cp_v,0,Wfd%nlmn_atm)
148    ABI_DT_MALLOC(Cp_c,(Wfd%natom,nspinor))
149    call pawcprj_alloc(Cp_c,0,Wfd%nlmn_atm)
150  end if
151 
152  if (inclvkb==1.and.usepaw==0) then
153    MSG_ERROR("inclvkb==1 not coded,using inclvkb==2")
154  end if
155  !
156  ! Calculate the matrix elements of ir in the IBZ.
157  ABI_MALLOC(ir_kibz,(3,lomo_min:max_band,lomo_min:max_band,Wfd%nkibz,nsppol))
158  ir_kibz=czero
159 
160  ABI_MALLOC(bbp_distrb, (Wfd%mband,Wfd%mband))
161 
162  do spin=1,nsppol
163    do ik_ibz=1,Wfd%nkibz
164     !
165     ! Distribute the (b,b') entries.
166     bbp_mask=.FALSE.; bbp_mask(lomo_spin(spin):max_band,lomo_spin(spin):max_band)=.TRUE.
167     call wfd_distribute_bbp(Wfd,ik_ibz,spin,"All",my_nbbp,bbp_distrb,bbp_mask=bbp_mask)
168     if (ALL(bbp_distrb/=my_rank)) CYCLE
169 
170     istwf_k = Wfd%istwfk(ik_ibz)
171     ABI_CHECK(istwf_k==1,"istwf_k/=1 not coded") ! KB stuff is missing.
172     npw_k = Wfd%npwarr(ik_ibz)
173     kg_k  => Wfd%Kdata(ik_ibz)%kg_k
174 
175     if (inclvkb/=0.and.usepaw==0) then
176       ! Prepare term i <n,k|[Vnl,r]|n"k>
177       call vkbr_init(vkbr,Cryst,Psps,inclvkb,istwf_k,npw_k,Kmesh%ibz(:,ik_ibz),kg_k)
178     end if
179 
180     ! Note: spinorial case is not coded therefore we work with ihrc(:,1).
181     ! TODO: The lower triangle can be Reconstructed by symmetry.
182     do ib_v=lomo_spin(spin),max_band ! Loop over bands
183       if ( ALL(bbp_distrb(ib_v,:)/=my_rank) ) CYCLE
184       ug_v => Wfd%Wave(ib_v,ik_ibz,spin)%ug
185       if (usepaw==1) then
186         call wfd_get_cprj(Wfd,ib_v,ik_ibz,spin,Cryst,Cp_v,sorted=.FALSE.)
187       end if
188 
189       do ib_c=lomo_spin(spin),max_band
190        if (bbp_distrb(ib_v,ib_c)/=my_rank) CYCLE
191        ug_c => Wfd%Wave(ib_c,ik_ibz,spin)%ug
192 
193        if (usepaw==0) then
194          ! Calculate matrix elements of i[H,r] for NC pseudopotentials.
195          ihrc = nc_ihr_comm(vkbr,cryst,psps,npw_k,nspinor,istwf_k,inclvkb,Kmesh%ibz(:,ik_ibz),ug_c,ug_v,kg_k)
196 
197        else
198          ! Matrix elements of i[H,r] for PAW.
199          call wfd_get_cprj(Wfd,ib_c,ik_ibz,spin,Cryst,Cp_c,sorted=.FALSE.)
200 
201          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)
202        end if
203        !
204        ! Save matrix elements of i*r in the IBZ
205        ediff = KS_Bst%eig(ib_c,ik_ibz,spin) - KS_BSt%eig(ib_v,ik_ibz,spin)
206        if (ABS(ediff)<tol16) ediff=tol6  ! Treat a possible degeneracy between v and c.
207        ir_kibz(:,ib_c,ib_v,ik_ibz,spin) = ihrc(:,1)/ediff
208 
209       end do !ib_c
210     end do !ib_v
211 
212     call vkbr_free(vkbr)
213    end do !spin
214  end do !ik_ibz
215 
216  ! Collect results on each node.
217  call xmpi_sum(ir_kibz,comm,ierr)
218 
219  ABI_FREE(bbp_distrb)
220 
221  if (usepaw==1) then
222    call pawcprj_free(Cp_v)
223    ABI_DT_FREE(Cp_v)
224    call pawcprj_free(Cp_c)
225    ABI_DT_FREE(Cp_c)
226  end if
227  !
228  ! ======================================================
229  ! ==== Calculate Fcv(kBZ) in the full Brilouin zone ====
230  ! ======================================================
231  !
232  ! Symmetrization of the matrix elements of the position operator.
233  ! <Sk b|r|Sk b'> = R <k b|r|k b'> + \tau \delta_{bb'}
234  !   where S is one of the symrec operations in reciprocal space, R is the
235  !   corresponding operation in real space, \tau being the associated fractional translations.
236  !
237  ! q.Mcv( Sk) =  S^{-1}q. Mcv(k)
238  ! q.Mcv(-Sk) = -S^{-1}q. CONJG(Mcv(k)) if time-reversal is used.
239 
240  b1=Cryst%gprimd(:,1)*two_pi
241  b2=Cryst%gprimd(:,2)*two_pi
242  b3=Cryst%gprimd(:,3)*two_pi
243 
244  opt_cvk = czero
245  do spin=1,nsppol
246    do ik_bz=1,nkbz
247     !
248     ! Get ik_ibz, and symmetries index from ik_bz.
249     call get_BZ_item(Kmesh,ik_bz,kbz,ik_ibz,isym_k,itim_k)
250 
251     mat_dp = DBLE(Cryst%symrec(:,:,isym_k))
252     call matrginv(mat_dp,3,3) ! Invert
253     qrot = (3-2*itim_k) * MATMUL(mat_dp,qpoint)
254 
255     do ib_v=lomo_spin(spin),max_band !  Loops over the bands C and V start
256       do ib_c=lomo_spin(spin),max_band
257         !if (ib_c==ib_v) CYCLE
258         emcvk = pdtqrc(qrot,ir_kibz(:,ib_c,ib_v,ik_ibz,spin),b1,b2,b3)
259         if (itim_k==2) emcvk = CONJG(emcvk)
260         opt_cvk(ib_c,ib_v,ik_bz,spin) = emcvk
261       end do !ib_c
262     end do !ib_v
263 
264    end do !ik_bz
265  end do !spin
266 
267  ABI_FREE(ir_kibz)
268 
269  call xmpi_barrier(comm)
270 
271 contains

ABINIT/m_wfd_optic [ Modules ]

[ Top ] [ Modules ]

NAME

  m_wfd_optic

FUNCTION

  Functions to compute optical matrix elements using the wavefunction descriptor.

COPYRIGHT

  Copyright (C) 2008-2018 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 .

PARENTS

CHILDREN

SOURCE

21 #if defined HAVE_CONFIG_H
22 #include "config.h"
23 #endif
24 
25 #include "abi_common.h"
26 
27 module m_wfd_optic
28 
29  use defs_basis
30  use m_errors
31  use m_abicore
32  use m_xmpi
33 
34  use defs_datatypes,      only : ebands_t, pseudopotential_type
35  use m_hide_lapack,       only : matrginv
36  use m_bz_mesh,           only : kmesh_t, get_BZ_item
37  use m_crystal,           only : crystal_t
38  use m_vkbr,              only : vkbr_t, vkbr_free, vkbr_init, nc_ihr_comm
39  use m_wfd,               only : wfd_t, wfd_get_cprj, wfd_distribute_bbp
40  use m_pawtab,            only : pawtab_type
41  use m_pawcprj,           only : pawcprj_type, pawcprj_alloc, pawcprj_free
42  use m_paw_hr,            only : pawhur_t, paw_ihr
43 
44  implicit none
45 
46  private

ABINIT/pdtqrc [ Functions ]

[ Top ] [ 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

PARENTS

SOURCE

289 pure function pdtqrc(R,C,b1,b2,b3)
290 
291 
292 !This section has been created automatically by the script Abilint (TD).
293 !Do not modify the following lines by hand.
294 #undef ABI_FUNC
295 #define ABI_FUNC 'pdtqrc'
296 !End of the abilint section
297 
298  implicit none
299 
300 !Arguments ------------------------------------
301 !arrays
302  real(dp),intent(in) :: R(3),b1(3),b2(3),b3(3)
303  complex(dpc),intent(in) :: C(3)
304  complex(dpc) :: pdtqrc
305 
306 !Local variables ------------------------------
307 !scalars
308  integer :: ii
309 
310 !************************************************************************
311 
312  pdtqrc=czero
313  do ii=1,3
314    pdtqrc = pdtqrc + (R(1)*b1(ii)+R(2)*b2(ii)+R(3)*b3(ii)) * &
315 &                    (C(1)*b1(ii)+C(2)*b2(ii)+C(3)*b3(ii))
316  end do
317 
318 end function pdtqrc