TABLE OF CONTENTS


ABINIT/calc_ucrpa [ Functions ]

[ Top ] [ Functions ]

NAME

 calc_ucrpa

FUNCTION

 Calculate the effective interaction in the correlated orbitals

COPYRIGHT

 Copyright (C) 1999-2018 ABINIT group (TApplencourt,BA)
 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.txt .

INPUTS

 npwe : number of plane wave for the dielectric constant
 npw : number of plane wave
 nomega  : number of frequencis
 bandinf,bandsup : kohn sham band
 optimisation : string for the optimisation
 Wfd:: MPI communicator
 mesh <kmesh_t>
    %nbz=Number of points in the BZ
    %nibz=Number of points in IBZ
    %kibz(,nibz)=k-point coordinates, irreducible Brillouin zone
    %kbz(3,nbz)=k-point coordinates, full Brillouin zone
    %ktab(nbz)= table giving for each k-point in the BZ (kBZ), the corresponding
    %ktabi(nbz)= for each k-point in the BZ defines whether inversion has to be considered
    %ktabp(nbz)= phase factor associated to tnons
 M1_q_m(bandinf:bandsup,bandinf:bandsup,npw,Qmesh%nibz): Oscillator strengh in Wannier basis
 rhot1_q_m(bandinf:bandsup,bandinf:bandsup,npw,Qmesh%nibz): Oscillator strengh
                                          multiplied by coulomb potential in Wannier basis

OUTPUT

NOTES

PARENTS

      sigma

CHILDREN

      affichage,checkk,cpu_time,get_bz_item,read_screening,wrtout
      xmpi_barrier,xmpi_sum

SOURCE

 88  subroutine calc_ucrpa(itypatcor,cryst,Kmesh,lpawu,M1_q_m,Qmesh,npwe,&
 89 & npw,nsym,rhot1_q_m,nomega,omegamin,omegamax,bandinf,bandsup,optimisation,ucvol,Wfd,fname)
 90 
 91  use defs_basis
 92  use m_abicore
 93  use m_xmpi
 94  use m_errors
 95 
 96  use m_io_tools,      only : open_file
 97  use m_wfd,           only : wfd_t
 98  use m_io_screening,  only : read_screening, em1_ncname
 99  use m_bz_mesh,       only : kmesh_t, get_BZ_item
100  use m_crystal,       only : crystal_t
101 
102 !This section has been created automatically by the script Abilint (TD).
103 !Do not modify the following lines by hand.
104 #undef ABI_FUNC
105 #define ABI_FUNC 'calc_ucrpa'
106 !End of the abilint section
107 
108  implicit none
109 !   _____            _
110 !  |_   _|          | |
111 !    | |  _ __  _ __  _   _| |_
112 !    | | | '_ \| '_ \| | | | __|
113 !   _| |_| | | | |_) | |_| | |_
114 !  |_____|_| |_| .__/ \__,_|\__|
115 !              | |
116 !              |_|
117 
118 !Arguments ------------------------------------
119  integer, intent(in)   :: itypatcor,lpawu,npw,npwe,nsym
120  integer, intent(in)   :: nomega
121  integer, intent(in)   :: bandinf
122  integer, intent(in)   :: bandsup
123  character(len=fnlen), intent(in) :: fname
124  character(len=*), intent(in) :: optimisation
125  real(dp), intent(in) :: ucvol,omegamin,omegamax
126 
127  type(wfd_t),intent(inout) :: Wfd
128  type(kmesh_t),intent(in) :: Kmesh,Qmesh
129  type(crystal_t),intent(in) :: Cryst
130  complex(dpc), intent(in) :: rhot1_q_m(cryst%nattyp(itypatcor),Wfd%nspinor,Wfd%nspinor,2*lpawu+1,2*lpawu+1,npw,Qmesh%nibz)
131  complex(dpc), intent(in) :: M1_q_m(cryst%nattyp(itypatcor),Wfd%nspinor,Wfd%nspinor,2*lpawu+1,2*lpawu+1,npw,Qmesh%nibz)
132 
133 !Local variables ------------------------------
134 !scalars
135  real(dp) :: x
136  real(dp) :: t1,t2
137  real(dp):: tol,eVnorme
138  complex(dpc) :: uu,jj
139 
140  complex :: nC,ualter
141 
142  integer :: iat,im_paral,iqalloc,ib1,ib2,m1,m2,m3,m4
143  integer :: ierr,ik_bz,ik_ibz,iq_ibz,i,iG1,iG2,iG,iiG,iomega,iomega1,ispinor1,ispinor2,ispinor3,ispinor4
144  integer :: lpawu_read,nkibz,nbband,nkbz,nprocs,nqalloc,nqibz,ms1,ms2,ms3,ms4,mbband,nspinor
145  integer :: isym_kgw,iik,unt
146  complex(dpc) ::ph_mkt
147 
148  logical :: wannier=.TRUE.
149  logical :: verbose=.FALSE.
150  logical :: bug=.FALSE.
151  logical :: lscr_one
152 
153  character(len=500) :: message
154 
155 !arrays
156  complex(dpc), allocatable :: V_m(:,:,:,:)
157  complex(dpc), allocatable :: U_m(:,:,:,:)
158 
159 ! complex(dpc), allocatable :: coeffW_BZ(:,:,:),coeffW_IBZ(:,:,:)
160  complex(dpc), allocatable :: rhot_q_m(:,:,:,:,:,:)
161  complex(dpc), allocatable :: rhot_q_m_npwe(:,:,:,:,:,:),trrho(:,:),sumrhorhoeps(:)
162  complex(gwpc), allocatable :: scr(:,:,:,:)
163 
164  real(dp),allocatable :: k_coord(:,:)!,k_coordIBZ(:,:)
165  real(dp),allocatable :: q_coord(:,:)
166  real(dp),allocatable:: normG(:)
167  complex(dpc),allocatable:: uomega(:),jomega(:)
168  real(dp),allocatable:: omega(:)
169 
170  integer,allocatable:: ikmq_bz_t(:,:)
171 
172  logical,allocatable :: bijection(:)
173 !************************************************************************
174 
175  write(message,*) ch10, '==== Calculation of the screened interaction ===='
176  call wrtout(std_out,message,'COLL')
177  call wrtout(ab_out,message,'COLL')
178  write(message,*) ""
179  call wrtout(std_out,message,'COLL')
180  call wrtout(ab_out,message,'COLL')
181  nkbz = Kmesh%nbz
182  nqibz= Qmesh%nibz
183  nspinor=Wfd%nspinor
184 
185  nbband=1+bandsup-bandinf
186 !  _  __            ____
187 ! | |/ /   ___     / __ \
188 ! | ' /   ( _ )   | |  | |
189 ! |  <    / _ \/\ | |  | |
190 ! | . \  | (_>  < | |__| |
191 ! |_|\_\  \___/\/  \___\_\
192 
193  write(message,*) "Read K and Q mesh"
194  call wrtout(std_out,message,'COLL')
195  call wrtout(ab_out,message,'COLL')
196 
197  ABI_ALLOCATE(k_coord,(nkbz,3))
198  ABI_ALLOCATE(q_coord,(nqibz,4))
199 
200 
201 !==Read k and q==!
202 !open(unit=2012,file='ikbz_COORD',form='formatted',status='unknown')
203 !read(2012,*) (ik_bz,k_coord(ik_bz,:),i=1,nkbz)
204 !close(2012)
205 
206  do ik_bz=1,nkbz
207    call get_BZ_item(Kmesh,ik_bz,k_coord(ik_bz,:),ik_ibz,isym_kgw,iik,ph_mkt)
208  end do
209 
210 ! open(unit=2012,file='iqbz_COORD',form='formatted',status='unknown')
211 ! read(2012,*)
212  do i=1,nqibz
213 !   read(2012,*) iq_ibz,q_coord(iq_ibz,:)
214    q_coord(i,1)=Qmesh%ibz(1,i)
215    q_coord(i,2)=Qmesh%ibz(2,i)
216    q_coord(i,3)=Qmesh%ibz(3,i)
217    q_coord(i,4)=Qmesh%wt(i)
218 !   if (iq_ibz > nqibz) then
219 !     write(message,*) iq_ibz,nqibz," Error on line",i,"Are you in iBZ ?"
220 !     call wrtout(std_out,message,'COLL')
221 !   end if
222  end do
223 ! close(2012)
224 
225 !==Bijection and array for k-q==!
226  ABI_ALLOCATE(bijection,(nkbz))
227  ABI_ALLOCATE(ikmq_bz_t,(nkbz,nqibz))
228  bijection(:)=.FALSE.
229  do ik_bz=1,nkbz
230    do iq_ibz=1,nqibz
231      ikmq_bz_t(ik_bz,iq_ibz)=findkmq(ik_bz,k_coord,q_coord(iq_ibz,:),nkbz)
232      if (ikmq_bz_t(ik_bz,iq_ibz)>nkbz.and.nsym==1) then
233        BUG=.TRUE.
234        write(message,*) "No K-Q for K/Q =",ik_bz,iq_ibz
235        MSG_ERROR(message)
236      end if
237      bijection(ikmq_bz_t(ik_bz,iq_ibz))=.TRUE.
238    end do
239 
240    if (count(bijection).NE.nqibz.and.nsym==1) then
241    BUG=.TRUE.
242    write(message,*) 'No bijection ',ik_bz
243    MSG_ERROR(message)
244    end if
245 
246    bijection(:)=.FALSE.
247  end do
248 
249  if (.NOT.BUG.and.nsym==1) then
250    write(message,*)  "Bijection Ok."
251    call wrtout(std_out,message,'COLL')
252  end if
253 !                                           _____
254 !                                          / ____|
255 !   _ __   ___  _ __ _ __ ___   ___       | |  __
256 !  | '_ \ / _ \| '__| '_ ` _ \ / _ \      | | |_ |
257 !  | | | | (_) | |  | | | | | |  __/      | |__| |
258 !  |_| |_|\___/|_|  |_| |_| |_|\___|       \_____|
259 
260  ABI_ALLOCATE(normG,(npw))
261  if (verbose) then
262    write(message,*) 'Read the potential and G norm'
263    call wrtout(std_out,message,'COLL')
264    if (open_file('normeG',message,newunit=unt,form='formatted',status='unknown') /= 0) then
265      MSG_ERROR(message)
266    end if
267    read(unt,*) (iiG,x,normG(iiG),iG=1,npw)
268    close(unt)
269    !!False norme for G=0 idd G is the inverse of the potential inverse du potentiel (q=0)
270    normG(1)=0
271  end if
272 
273 !========================================================================
274 !------------------------------------------------------------------------
275 
276 !  FIRST PART OF THE ROUTINE: USE M_G^(nn')(q,k) to do formal checks.
277 !                             USE rhot_q_n to do compute bare interaction
278 
279 !------------------------------------------------------------------------
280 !========================================================================
281 
282 !                               _
283 !                              ( )
284 !   _ __ ___        _ __  _ __ |/
285 !  | '_ ` _ \      | '_ \| '_ \
286 !  | | | | | |     | | | | | | |
287 !  |_| |_| |_|     |_| |_|_| |_|
288 
289 !==========================================================
290 !==========================================================
291 ! tol=1E-1
292 ! tolerance for the normalization of wfc: should be around 0.01.
293  tol = 1 ! very large for test.
294  write(message,*) 'Check the norm of M'
295  call wrtout(std_out,message,'COLL')
296  write(message,*) 'Tolerance :',tol
297  call wrtout(std_out,message,'COLL')
298 !  __      __
299 !  \ \    / /
300 !   \ \  / /      _ __
301 !    \ \/ /      | '_ \
302 !     \  /       | | | |
303 !      \/        |_| |_|
304 !==========================================================================
305 !==Compute V_{n,n'}: bare interaction in the KS basis <- rhot_q_n -> V_n
306 !==========================================================================
307 !==========================================================
308 !==Compute V_{n,n'}
309 !==========================================================
310  if(verbose) then
311    write(message,*) ""
312    call wrtout(std_out,message,'COLL')
313    call wrtout(ab_out,message,'COLL')
314    write(message,*)  "==Calcul of the bare kohn-sham interaction V n=="
315    call wrtout(std_out,message,'COLL')
316  endif
317  tol=1E+1
318 
319  if (.NOT.wannier) RETURN
320 
321 !========================================================================
322 !------------------------------------------------------------------------
323 
324 !  SECOND PART OF THE ROUTINE: Read Wannier projections
325 
326 !------------------------------------------------------------------------
327 !========================================================================
328 !
329 ! \ \        / /       (_)
330 !  \ \  /\  / /_ _ _ __  _ __  _  ___ _ __
331 !   \ \/  \/ / _` | '_ \| '_ \| |/ _ \ '__|
332 !    \  /\  / (_| | | | | | | | |  __/ |
333 !     \/  \/ \__,_|_| |_|_| |_|_|\___|_|
334 
335 !==========================================================
336 !==========================================================
337 !== Read Wannier coefficient in forlb.ovlp
338 !==========================================================
339  write(message,*) ""
340  call wrtout(std_out,message,'COLL')
341  call wrtout(ab_out,message,'COLL')
342  write(message,*) "Read wannier in iBZ"
343  call wrtout(ab_out,message,'COLL')
344  call wrtout(std_out,message,'COLL')
345  nkibz=Kmesh%nibz
346 
347 !Read "l"
348  if (open_file('forlb.ovlp',message,newunit=unt,form='formatted',status='unknown') /= 0) then
349    MSG_ERROR(message)
350  end if
351  rewind(unt)
352  read(unt,*) message, lpawu_read
353  read(unt,*) message, ib1, ib2
354  close(unt)
355  if(ib1/=bandinf.and.ib2/=bandsup) then
356    write(message,*) "Error with bands",ib1,bandinf,ib2,bandsup
357    MSG_ERROR(message)
358  endif
359 !!Read the bandinf, bandinf redondance information
360  mbband=2*lpawu_read+1
361 
362 !*******************************************************
363 !if (3==4) then
364 !!  USELESS START

ABINIT/m_calc_ucrpa [ Modules ]

[ Top ] [ Modules ]

NAME

  m_calc_ucrpa

FUNCTION

 Calculate the effective interaction in the correlated orbitals

COPYRIGHT

 Copyright (C) 2006-2018 ABINIT group (BAmadon)
 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

OUTPUT

PARENTS

CHILDREN

SOURCE

24 #if defined HAVE_CONFIG_H
25 #include "config.h"
26 #endif
27 
28 
29 #include "abi_common.h"
30 
31 MODULE m_calc_ucrpa
32 
33  use defs_basis
34  implicit none
35 
36  private
37 
38  public :: calc_ucrpa