TABLE OF CONTENTS


ABINIT/m_oper [ Modules ]

[ Top ] [ Modules ]

NAME

  m_oper

FUNCTION

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_oper
32 
33  use defs_basis
34  use m_abicore
35  use m_errors
36 
37  use m_matlu,    only : matlu_type
38  use m_hide_lapack,  only : xginv
39 
40  implicit none
41 
42  private
43 
44  public :: init_oper
45  public :: diff_oper
46  public :: destroy_oper
47  public :: print_oper
48  public :: inverse_oper
49  public :: loc_oper
50  public :: identity_oper
51  public :: copy_oper
52  public :: trace_oper
53  public :: upfold_oper
54  public :: prod_oper

m_oper/copy_oper [ Functions ]

[ Top ] [ m_oper ] [ Functions ]

NAME

 copy_oper

FUNCTION

INPUTS

OUTPUT

PARENTS

      hybridization_asymptotic_coefficient,m_green

CHILDREN

      prod_matlu

SOURCE

302 subroutine copy_oper(oper1,oper2)
303 
304  use defs_basis
305  use m_matlu, only : copy_matlu
306  use m_errors
307 
308 !This section has been created automatically by the script Abilint (TD).
309 !Do not modify the following lines by hand.
310 #undef ABI_FUNC
311 #define ABI_FUNC 'copy_oper'
312 !End of the abilint section
313 
314  implicit none
315 
316 !Arguments ------------------------------------
317 !type
318  type(oper_type),intent(in) :: oper1
319  type(oper_type),intent(inout) :: oper2 !vz_i
320 
321 !oper variables-------------------------------
322  integer ::  ib, ib1, ikpt, isppol
323 ! *********************************************************************
324  DBG_ENTER("COLL")
325  if(oper1%has_opermatlu==1.and.oper2%has_opermatlu==1)  then
326    call copy_matlu(oper1%matlu,oper2%matlu,oper1%natom)
327  endif
328 
329  if(oper1%has_operks==1.and.oper2%has_operks==1) then
330    do isppol=1,oper1%nsppol
331      do ikpt=1,oper1%nkpt
332        do ib=1,oper1%mbandc
333          do ib1=1,oper1%mbandc
334            oper2%ks(isppol,ikpt,ib,ib1)=oper1%ks(isppol,ikpt,ib,ib1)
335          enddo
336        enddo
337      enddo
338    enddo
339  endif
340 
341  DBG_EXIT("COLL")
342 end subroutine copy_oper

m_oper/destroy_oper [ Functions ]

[ Top ] [ m_oper ] [ Functions ]

NAME

 destroy_oper

FUNCTION

  deallocate oper

INPUTS

  oper

OUTPUT

PARENTS

      hubbard_one,hybridization_asymptotic_coefficient,m_green,m_self
      outscfcv,psichi_renormalization,qmc_prep_ctqmc,vtorho

CHILDREN

      prod_matlu

SOURCE

240 subroutine destroy_oper(oper)
241 
242  use defs_basis
243  use m_crystal, only : crystal_t
244  use m_matlu, only : destroy_matlu
245  use m_errors
246 
247 !This section has been created automatically by the script Abilint (TD).
248 !Do not modify the following lines by hand.
249 #undef ABI_FUNC
250 #define ABI_FUNC 'destroy_oper'
251 !End of the abilint section
252 
253  implicit none
254 
255 !Arguments ------------------------------------
256 !scalars
257  type(oper_type),intent(inout) :: oper
258 !local variables-------------------------------
259  character(len=500) :: message
260 
261 !! *********************************************************************
262  DBG_ENTER("COLL")
263  if(oper%has_opermatlu==1) then
264    call destroy_matlu(oper%matlu,oper%natom)
265  else
266    message = " Operator is not defined to be used in destroy_oper"
267    MSG_ERROR(message)
268  endif
269  if ( allocated(oper%matlu))  then
270    ABI_DATATYPE_DEALLOCATE(oper%matlu)
271    oper%has_opermatlu=0
272  endif
273  if ( allocated(oper%ks)) then
274    ABI_DEALLOCATE(oper%ks)
275    oper%has_operks=0
276  endif
277  oper%wtk => null()
278 !  no deallocation for wtk: wtk is an explicit pointer
279 
280  DBG_EXIT("COLL")
281 end subroutine destroy_oper

m_oper/diff_oper [ Functions ]

[ Top ] [ m_oper ] [ Functions ]

NAME

 diff_oper

FUNCTION

 Compute a norm of the differences between two occupations matrices.

INPUTS

  cryst_struc <type(crystal_t)>=crystal structure data
  occup1 <type(oper_type)>= occupations
  occup2 <type(oper_type)>= occupations
  option : option for printing (if 1 assume data are related to lda only)
  paw_dmft  <type(paw_dmft_type)>= paw+dmft related data

OUTPUT

PARENTS

      m_dmft,psichi_renormalization

CHILDREN

      prod_matlu

SOURCE

 956 subroutine diff_oper(char1,char2,occup1,occup2,option,toldiff)
 957 
 958  use defs_basis
 959  use m_paw_dmft, only : paw_dmft_type
 960  use m_crystal, only : crystal_t
 961  use m_matlu, only : diff_matlu
 962  use m_errors
 963 
 964 !This section has been created automatically by the script Abilint (TD).
 965 !Do not modify the following lines by hand.
 966 #undef ABI_FUNC
 967 #define ABI_FUNC 'diff_oper'
 968 !End of the abilint section
 969 
 970  implicit none
 971 
 972 !Arguments ------------------------------------
 973 !type
 974  type(oper_type), intent(in) :: occup1,occup2
 975  integer, intent(in) :: option
 976  real(dp), intent(in) :: toldiff
 977  character(len=*), intent(in) :: char1,char2
 978 !local variables-------------------------------
 979  integer :: mbandc,nkpt
 980  character(len=500) :: message
 981 ! *********************************************************************
 982 
 983  DBG_ENTER("COLL")
 984  if(occup1%has_opermatlu==0.or.occup2%has_opermatlu==0) then
 985    message = " Operators are not defined to be used in diff_oper"
 986    MSG_ERROR(message)
 987  endif
 988  mbandc   = occup1%mbandc
 989  nkpt    = occup1%nkpt
 990  if(occup1%nkpt/=occup2%nkpt) then
 991   write(message,'(a,2x,2i9)')' Operators are not equals',occup1%nkpt,occup2%nkpt
 992   MSG_ERROR(message)
 993  endif
 994 
 995  call diff_matlu(char1,char2,occup1%matlu,&
 996 & occup2%matlu,occup1%natom,option,toldiff)
 997 ! if(option==1) then
 998 !  toldiff=tol4
 999 !  if( matludiff < toldiff ) then
1000 !   write(message,'(6a,e12.4,a,e12.4)') ch10,&
1001 !&   '   Differences between ',trim(char1),' and ',trim(char2),' is small enough:',&
1002 !&   matludiff,'is lower than',toldiff
1003 !   call wrtout(std_out,message,'COLL')
1004 !  else
1005 !   write(message,'(6a,e12.4,a,e12.4)') ch10,&
1006 !&   '   Error: Differences between ',trim(char1),' and ',trim(char2),' is too large:',&
1007 !&   matludiff,'is large than',toldiff
1008 !   call wrtout(std_out,message,'COLL')
1009 !   call abi_abort('COLL')
1010 !  endif
1011 ! endif
1012 ! call abi_abort('COLL')
1013 
1014  DBG_EXIT("COLL")
1015 end subroutine diff_oper

m_oper/identity_oper [ Functions ]

[ Top ] [ m_oper ] [ Functions ]

NAME

 identity_oper

FUNCTION

INPUTS

OUTPUT

PARENTS

      psichi_renormalization

CHILDREN

      prod_matlu

SOURCE

862 subroutine identity_oper(oper,option)
863 
864  use defs_basis
865  use m_crystal, only : crystal_t
866  use m_paw_dmft, only : paw_dmft_type
867  use m_errors
868 
869 !This section has been created automatically by the script Abilint (TD).
870 !Do not modify the following lines by hand.
871 #undef ABI_FUNC
872 #define ABI_FUNC 'identity_oper'
873 !End of the abilint section
874 
875  implicit none
876 
877 !Arguments ------------------------------------
878 !type
879  integer, intent(in):: option
880  type(oper_type),intent(inout) :: oper
881 !oper variables-------------------------------
882  integer :: iatom,ib,ikpt,ispinor,isppol,im
883  integer :: natom,mbandc,ndim,nkpt,nspinor,nsppol
884  character(len=500) :: message
885 ! *********************************************************************
886 
887  DBG_ENTER("COLL")
888 
889  if(((option==1.or.option==3).and.(oper%has_opermatlu==0)).or.&
890 &   ((option==2.or.option==3).and.(oper%has_operks==0))) then
891    message = " Options in identity_oper are not coherent with definitions of this operator"
892    MSG_ERROR(message)
893  endif
894  nkpt=oper%nkpt
895  nsppol=oper%nsppol
896  mbandc=oper%mbandc
897  natom=oper%natom
898  nspinor=oper%nspinor
899  oper%ks=czero
900  do iatom=1,natom
901   oper%matlu(iatom)%mat= czero
902  enddo
903 
904  if(option==1.or.option==3) then
905    do isppol=1,nsppol
906      do ikpt=1,nkpt
907        do ib=1,mbandc
908          oper%ks(isppol,ikpt,ib,ib)=cone
909        enddo ! ib
910      enddo ! ikpt
911    enddo ! isppol
912 
913  else if (option==2.or.option==3) then
914    do iatom=1,natom
915      if(oper%matlu(iatom)%lpawu.ne.-1) then
916        ndim=2*oper%matlu(iatom)%lpawu+1
917        do isppol=1,nsppol
918          do im=1,ndim
919            do ispinor=1,nspinor
920              oper%matlu(iatom)%mat(im,im,isppol,ispinor,ispinor)= cone
921            enddo ! ispinor
922          enddo ! im
923        enddo
924      endif ! lpawu
925    enddo ! iatom
926  endif
927 
928  DBG_EXIT("COLL")
929 end subroutine identity_oper

m_oper/init_oper [ Functions ]

[ Top ] [ m_oper ] [ Functions ]

NAME

 init_oper

FUNCTION

  Allocate variables used in type oper_type.

INPUTS

 OUTPUTS
 oper  = operator of type oper_type

PARENTS

      hubbard_one,hybridization_asymptotic_coefficient,m_green,m_self
      outscfcv,psichi_renormalization,qmc_prep_ctqmc,vtorho

CHILDREN

      prod_matlu

SOURCE

132 subroutine init_oper(paw_dmft,oper,nkpt,wtk,opt_ksloc)
133 
134  use defs_basis
135  use m_matlu, only : init_matlu
136  use m_paw_dmft, only : paw_dmft_type
137  use m_errors
138 
139 !This section has been created automatically by the script Abilint (TD).
140 !Do not modify the following lines by hand.
141 #undef ABI_FUNC
142 #define ABI_FUNC 'init_oper'
143 !End of the abilint section
144 
145  implicit none
146 
147 !Arguments ------------------------------------
148 !scalars
149  integer, optional, intent(in) :: nkpt
150  integer, optional, intent(in) :: opt_ksloc
151 !type
152  type(paw_dmft_type), intent(in) :: paw_dmft
153  type(oper_type), intent(inout) :: oper
154 !arrays
155  real(dp), pointer, optional :: wtk(:)
156 !oper variables ------------------------------------
157  integer :: iatom,optksloc
158 
159 !************************************************************************
160  DBG_ENTER("COLL")
161  if(present(opt_ksloc)) then
162    optksloc=opt_ksloc
163  else
164    optksloc=3
165  endif
166 
167  if(optksloc/=3) then
168     ! FIXME: empty line!
169  endif
170 
171  oper%has_operks=0
172  oper%has_opermatlu=0
173 ! ===================
174 !  Integers
175 ! ===================
176  oper%nsppol=paw_dmft%nsppol
177  oper%nspinor=paw_dmft%nspinor
178  oper%mbandc=paw_dmft%mbandc
179  oper%natom=paw_dmft%natom
180 
181 ! ===================
182 !  KS variables
183 ! ===================
184  if(optksloc==1.or.optksloc==3) then
185    if(.not.present(nkpt)) then
186      oper%nkpt=paw_dmft%nkpt
187    else
188      oper%nkpt=nkpt
189    endif
190 ! allocate(oper%wtk(oper%nkpt))
191    if(.not.present(wtk)) then
192      oper%wtk=>paw_dmft%wtk
193    else
194      oper%wtk=>wtk
195    endif
196    ABI_ALLOCATE(oper%ks,(paw_dmft%nsppol,oper%nkpt,paw_dmft%mbandc,paw_dmft%mbandc))
197    oper%has_operks=1
198    oper%ks=czero
199  endif
200 
201 ! ===================
202 !  matlu variables
203 ! ===================
204  if(optksloc==2.or.optksloc==3) then
205    oper%has_opermatlu=0
206    ABI_DATATYPE_ALLOCATE(oper%matlu,(oper%natom))
207    oper%has_opermatlu=1
208    call init_matlu(oper%natom,paw_dmft%nspinor,paw_dmft%nsppol,paw_dmft%lpawu,oper%matlu)
209    do iatom=1,oper%natom
210     oper%matlu(iatom)%mat=czero
211    enddo
212  endif
213 
214  DBG_EXIT("COLL")
215 
216 end subroutine init_oper

m_oper/inverse_oper [ Functions ]

[ Top ] [ m_oper ] [ Functions ]

NAME

 inverse_oper

FUNCTION

  Compute the inverse of the operator either in the KS space or in the
  correlated subspace.

INPUTS

  oper <type(oper_type)>= operator
  paw_dmft <type(paw_dmft_type)>
  option= integer

OUTPUT

  oper <type(oper_type)>= operator inverted

PARENTS

      dyson,m_green,qmc_prep_ctqmc

CHILDREN

      prod_matlu

SOURCE

524 subroutine inverse_oper(oper,option,prtopt,procb,iproc)
525 
526  use defs_basis
527  use m_crystal, only : crystal_t
528  use m_paw_dmft, only : paw_dmft_type
529  use m_matlu, only : inverse_matlu
530  use m_errors
531 
532 !This section has been created automatically by the script Abilint (TD).
533 !Do not modify the following lines by hand.
534 #undef ABI_FUNC
535 #define ABI_FUNC 'inverse_oper'
536 !End of the abilint section
537 
538  implicit none
539 
540 !Arguments ------------------------------------
541 !type
542  integer, intent(in):: option
543  integer, intent(in) :: prtopt
544  type(oper_type),intent(inout) :: oper
545  integer, optional, intent(in) ::  iproc
546  integer, optional, intent(in) :: procb(oper%nkpt)
547 !oper variables-------------------------------
548  integer :: ikpt,isppol
549  complex(dpc), allocatable :: matrix(:,:)
550  character(len=500) :: message
551  integer :: paral
552  integer, allocatable :: procb2(:)
553 !todo_ba: prb with gwpc here: necessary for matcginv but should be dpc
554 ! *********************************************************************
555  DBG_ENTER("COLL")
556  if(option==2.or.option==3) then
557    ABI_ALLOCATE(procb2,(oper%nkpt))
558  endif
559 !  if option=1 do inversion in local space
560 !  if option=2 do inversion in KS band space
561 !  if option=3 do both
562  if(present(procb).and.present(iproc)) then
563    paral=1
564    procb2=procb
565  else
566    paral=0
567  endif
568 
569  if(((option==1.or.option==3).and.(oper%has_opermatlu==0)).or.&
570 &   ((option==2.or.option==3).and.(oper%has_operks==0))) then
571    message = " Options are not coherent with definitions of this operator"
572    MSG_ERROR(message)
573  endif
574 
575  if(option==1.or.option==3) then
576    call inverse_matlu(oper%matlu,oper%natom,prtopt)
577  else if(option==2.or.option==3) then
578    ABI_ALLOCATE(matrix,(oper%mbandc,oper%mbandc))
579      do isppol=1,oper%nsppol
580        do ikpt=1,oper%nkpt
581         if ((paral==1.and.(procb2(ikpt)==iproc)).or.(paral==0)) then
582           matrix(:,:)=oper%ks(isppol,ikpt,:,:)
583 !          write(std_out,*) "isppol,ikpt",isppol,ikpt,m
584 !          write(std_out,*) "isppol,ikpt",matrix
585          !call matcginv_dpc(matrix,oper%mbandc,oper%mbandc)
586          call xginv(matrix,oper%mbandc)
587          oper%ks(isppol,ikpt,:,:)=matrix(:,:)
588         endif
589        enddo ! ikpt
590      enddo ! isppol
591    ABI_DEALLOCATE(matrix)
592  endif
593 
594  if(option==2.or.option==3) then
595    ABI_DEALLOCATE(procb2)
596  endif
597  DBG_EXIT("COLL")
598 end subroutine inverse_oper

m_oper/loc_oper [ Functions ]

[ Top ] [ m_oper ] [ Functions ]

NAME

 loc_oper

FUNCTION

INPUTS

OUTPUT

PARENTS

      compute_levels,hybridization_asymptotic_coefficient,m_green
      psichi_renormalization

CHILDREN

      prod_matlu

SOURCE

620 subroutine loc_oper(oper,paw_dmft,option,jkpt,procb,iproc)
621 
622  use defs_basis
623  use m_paw_dmft, only : paw_dmft_type
624  use m_errors
625 
626 !This section has been created automatically by the script Abilint (TD).
627 !Do not modify the following lines by hand.
628 #undef ABI_FUNC
629 #define ABI_FUNC 'loc_oper'
630 !End of the abilint section
631 
632  implicit none
633 
634 !Arguments ------------------------------------
635 !type
636  integer, intent(in):: option
637  integer, optional, intent(in):: jkpt
638  type(oper_type),intent(inout) :: oper
639  type(paw_dmft_type), intent(in) :: paw_dmft
640  integer, optional, intent(in) ::  iproc
641  integer, optional, intent(in) :: procb(oper%nkpt)
642 !oper variables-------------------------------
643  integer :: iatom,ib,ib1,ikpt,ikpt1,ispinor,ispinor1,isppol,im1,im
644  integer :: natom,mbandc,ndim,nkpt,nspinor,nsppol,paral
645  character(len=500) :: message
646  integer, allocatable :: procb2(:)
647  logical lvz  !vz_d
648 ! *********************************************************************
649  DBG_ENTER("COLL")
650  ABI_ALLOCATE(procb2,(oper%nkpt))
651  if((oper%has_opermatlu==0).or.(oper%has_operks==0)) then
652    message = " Operator is not defined to be used in loc_oper"
653    MSG_ERROR(message)
654  endif
655  if(present(procb).and.present(iproc)) then
656    paral=1
657    procb2=procb
658  else
659    paral=0
660  endif
661 
662  if(option<0) then
663  endif
664  nkpt=oper%nkpt
665  natom=oper%natom
666  nsppol=oper%nsppol
667  mbandc=oper%mbandc
668  nspinor=oper%nspinor
669 
670 
671  do iatom=1,natom
672    oper%matlu(iatom)%mat=czero
673  enddo
674  do isppol=1,nsppol
675   do ikpt=1,nkpt
676    ikpt1=ikpt
677    if(present(jkpt)) ikpt1=jkpt
678    lvz=paral==0  !vz_d
679    if(present(iproc)) lvz=lvz.or.(paral==1.and.(procb2(ikpt1)==iproc))  !vz_d
680 !  if ((paral==1.and.(procb2(ikpt1)==iproc)).or.(paral==0)) then    !vz_d
681    if(lvz) then !vz_d
682    do ib=1,mbandc
683     do ib1=1,mbandc
684      do iatom=1,natom
685       if(oper%matlu(iatom)%lpawu.ne.-1) then
686       ndim=2*oper%matlu(iatom)%lpawu+1
687        do im=1,ndim
688         do im1=1,ndim
689          do ispinor=1,nspinor
690           do ispinor1=1,nspinor
691             if (im1 == im .and. im1 == 1 .and. ib1 == ib .and. iatom == 1) then
692 
693             end if
694             oper%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)=     &
695 &            oper%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)+    &
696 &            paw_dmft%psichi(isppol,ikpt1,ib,ispinor,iatom,im)*        &
697 &            conjg(paw_dmft%psichi(isppol,ikpt1,ib1,ispinor1,iatom,im1))* &
698 !false&            paw_dmft%psichi(isppol,ikpt1,ib1,ispinor1,iatom,im1)*        &
699 !false&            conjg(paw_dmft%psichi(isppol,ikpt1,ib,ispinor,iatom,im))* &
700 &            oper%ks(isppol,ikpt,ib,ib1)*oper%wtk(ikpt)
701 ! one  could suppress wtk here if present(jkpt)
702 ! ks(ib,ib1)=ks(ib1,ib) -> ib and ib1 can be underchanged !
703           enddo ! ispinor1
704          enddo ! ispinor
705         enddo ! im1
706        enddo ! im
707       endif
708      enddo ! iatom
709     enddo ! ib
710    enddo ! ib
711    endif
712   enddo ! ikpt
713  enddo ! isppol
714  ABI_DEALLOCATE(procb2)
715 
716 
717 
718 
719  DBG_EXIT("COLL")
720 end subroutine loc_oper

m_oper/oper_type [ Types ]

[ Top ] [ m_oper ] [ Types ]

NAME

  oper_type

FUNCTION

  This structured datatype contains the necessary data

SOURCE

 66  type, public :: oper_type ! for each atom
 67 
 68 !  integer :: maxlpawu         ! Number of correlated atoms
 69 !
 70 !  integer :: mband
 71 !  ! Number of bands
 72 !
 73 !
 74   integer :: nkpt
 75   ! Number of k-point in the IBZ.
 76 !
 77   integer :: natom
 78 !
 79   integer :: mbandc
 80 !  ! Total number of bands in the Kohn-Sham Basis for PAW+DMFT
 81 !
 82   integer :: nspinor
 83 !
 84   integer :: nsppol
 85 
 86   integer :: has_operks
 87 
 88   integer :: has_opermatlu
 89 
 90   character(len=12) :: whichoper
 91   ! describe the type of operator computed (LDA, DMFT, KS..)
 92 
 93 !  ! Polarisation
 94   type(matlu_type), allocatable :: matlu(:)
 95 !   Local projection on correlated orbitals
 96 
 97   complex(dpc), allocatable :: ks(:,:,:,:)
 98 !   In the KS basis  (nsppol,nkpt,nband,nband)
 99 
100   real(dp), pointer :: wtk(:) => null()
101 
102  end type oper_type

m_oper/print_oper [ Functions ]

[ Top ] [ m_oper ] [ Functions ]

NAME

 print_oper

FUNCTION

INPUTS

 option= 1
         2
         3
         4
  below  5: write diagonal part of KS occupation matrix
         6
         7
  above  8: write all elements of KS occup. matrix.
         9

OUTPUT

PARENTS

      m_green,m_self

CHILDREN

      prod_matlu

SOURCE

373 subroutine print_oper(oper,option,paw_dmft,prtopt)
374 
375  use defs_basis
376  use m_matlu, only : print_matlu
377  use m_paw_dmft, only : paw_dmft_type
378  use m_errors
379 
380 !This section has been created automatically by the script Abilint (TD).
381 !Do not modify the following lines by hand.
382 #undef ABI_FUNC
383 #define ABI_FUNC 'print_oper'
384 !End of the abilint section
385 
386  implicit none
387 
388 !Arguments ------------------------------------
389 !type
390  type(paw_dmft_type), intent(in) :: paw_dmft
391  type(oper_type),intent(in) :: oper
392  integer, intent(in) :: option,prtopt
393 
394 !oper variables-------------------------------
395  integer :: ib,ib1,ikpt,isppol,nkptr,iband1,iband2
396  character(len=2000) :: message
397  logical :: ximag
398  real(dp) :: maximag(3)
399 ! *********************************************************************
400  DBG_ENTER("COLL")
401  maximag(:)=zero
402 
403  if(oper%has_opermatlu==1) then
404    write(message,'(2a)') ch10,'   = In the atomic basis'
405    call wrtout(std_out,message,'COLL')
406    call print_matlu(oper%matlu,oper%natom,prtopt)
407  endif
408  if(oper%has_operks==1) then
409    write(message,'(2a)') ch10,'   = In the KS basis'
410    call wrtout(std_out,message,'COLL')
411 
412 !todo_ba complete print_out
413    iband1=1
414    iband2=oper%mbandc
415 !   do ib=1,oper%mbandc
416 !     if(-(paw_dmft%eigen_lda(1,1,ib)+paw_dmft%fermie).ge.0.3) iband1=ib
417 !     if( (paw_dmft%eigen_lda(1,1,ib)-paw_dmft%fermie).le.0.3) iband2=ib
418 !   enddo
419 
420    ximag=.false.
421    if((abs(prtopt)>=3.and.((option<5).or.(option>8))).and.oper%has_operks==1) then
422 !     write(message,'(x,a,a,i4,2x,a)') ch10,'  -KS states'
423 !     call wrtout(std_out,message,'COLL')
424      do isppol=1,paw_dmft%nsppol
425        write(message, '(a,3x,a,i4)') ch10,"--isppol--",isppol
426        call wrtout(std_out,message,'COLL')
427        nkptr=oper%nkpt
428        nkptr=min(oper%nkpt,4)
429        write(message, '(2a)') ch10,&
430 &       "   - ( in the following only the value for the first k-points are printed)"
431        call wrtout(std_out,message,'COLL')
432        do ikpt=1,nkptr
433          if(option<5) then
434            write(message, '(2a,i4,2x,f14.5,a)') ch10,&
435 &           "   -k-pt--",ikpt,oper%wtk(ikpt),"(<-weight(k-pt))"
436            call wrtout(std_out,message,'COLL')
437          else if(abs(prtopt)>=4.or.option>8) then
438            write(message, '(2a,i5,a,i5,a,i5)') ch10,"  Writes occupations for k-pt",&
439 &           ikpt, "and between bands",&
440 &           iband1," and",iband2
441            call wrtout(std_out,message,'COLL')
442          endif
443          do ib=1,oper%mbandc
444            if(option<5) then
445              if(abs(aimag(oper%ks(isppol,ikpt,ib,ib))).ge.tol10) then
446                write(message, '(a,i4,e14.5,3x,e14.5,3x,e21.14)') "   -iband--",ib,&
447 &               paw_dmft%eigen_lda(isppol,ikpt,ib),oper%ks(isppol,ikpt,ib,ib)
448                call wrtout(std_out,message,'COLL')
449              else
450                write(message, '(a,i4,e14.5,3x,e14.5)') "   -iband--",ib,&
451 &               paw_dmft%eigen_lda(isppol,ikpt,ib),real(oper%ks(isppol,ikpt,ib,ib))
452                call wrtout(std_out,message,'COLL')
453              endif
454            endif
455            if(abs(prtopt)>=4.or.option>8.and.ib>=iband1.and.ib<=iband2) then
456              write(message, '(2000(f8.3))') &
457 &               (real(oper%ks(isppol,ikpt,ib,ib1)),ib1=iband1,iband2)
458              call wrtout(std_out,message,'COLL')
459              write(message, '(2000(f8.3))') &
460 &               (aimag(oper%ks(isppol,ikpt,ib,ib1)),ib1=iband1,iband2)
461              call wrtout(std_out,message,'COLL')
462              write(message, '(2000(f8.3))') &
463 &               (sqrt(real(oper%ks(isppol,ikpt,ib,ib1))**2+aimag(oper%ks(isppol,ikpt,ib,ib1))**2),ib1=iband1,iband2)
464              call wrtout(std_out,message,'COLL')
465 !   to write imaginary part
466 !             write(message, '(1000(2f9.3,2x))') &
467 !&               (real(oper%ks(isppol,ikpt,ib,ib1)),imag(oper%ks(isppol,ikpt,ib,ib1)),ib1=iband1,iband2)
468 !             call wrtout(std_out,message,'COLL')
469            endif ! prtopt>=20
470            do ib1=1,oper%mbandc
471              if(abs(aimag(oper%ks(isppol,ikpt,ib,ib1)))>max(tol10,maximag(1))) then
472                ximag=.true.
473                maximag(1)=aimag(oper%ks(isppol,ikpt,ib,ib1))
474              endif
475            enddo
476          enddo ! ib
477        enddo ! ikpt
478      enddo ! isppol
479    else
480     write(message,'(5x,a,i10,a)') '(not written)'
481     call wrtout(std_out,message,'COLL')
482    endif
483    if(ximag) then
484      write(message, '(3a,e12.4,a)')"Occupations are imaginary !",ch10, &
485 &    "  Maximal value is ", maximag(1), ch10
486      MSG_WARNING(message)
487    endif
488  else if(abs(prtopt)>=3.and.((option<5).or.(option>8))) then
489    write(message, '(2a)') ch10," Prb with options and has_operks in print_oper"
490    call wrtout(std_out,message,'COLL')
491  endif ! if oper%has_operks
492 ! write(message, '(2a)') ch10," end print_oper"
493 !     call wrtout(std_out,message,'COLL')
494 
495 
496  DBG_EXIT("COLL")
497 end subroutine print_oper

m_oper/prod_oper [ Functions ]

[ Top ] [ m_oper ] [ Functions ]

NAME

 prod_oper

FUNCTION

INPUTS

OUTPUT

PARENTS

      hybridization_asymptotic_coefficient

CHILDREN

      prod_matlu

SOURCE

1111 subroutine prod_oper(oper1,oper2,oper3,opt_ksloc)
1112 
1113  use defs_basis
1114  use m_errors
1115  use m_matlu, only : prod_matlu
1116 
1117 !This section has been created automatically by the script Abilint (TD).
1118 !Do not modify the following lines by hand.
1119 #undef ABI_FUNC
1120 #define ABI_FUNC 'prod_oper'
1121 !End of the abilint section
1122 
1123  implicit none
1124 
1125 !Arguments ------------------------------------
1126 !type
1127  type(oper_type),intent(in) :: oper1
1128  type(oper_type),intent(in) :: oper2
1129  type(oper_type),intent(inout) :: oper3
1130  integer :: opt_ksloc
1131 
1132 !oper variables-------------------------------
1133  integer ::  ib, ib1, ib2, ikpt, isppol
1134 ! *********************************************************************
1135  DBG_ENTER("COLL")
1136  if(opt_ksloc==2) then
1137    if(oper1%has_opermatlu==1.and.oper2%has_opermatlu==1.and.oper3%has_opermatlu==1)  then
1138      call prod_matlu(oper1%matlu,oper2%matlu,oper3%matlu,oper1%natom) !
1139    endif
1140  endif
1141 
1142  if(opt_ksloc==1) then
1143    if(oper1%has_operks==1.and.oper2%has_operks==1.and.oper3%has_operks==1) then
1144      do isppol=1,oper1%nsppol
1145        do ikpt=1,oper1%nkpt
1146          do ib=1,oper1%mbandc
1147            do ib1=1,oper1%mbandc
1148              oper3%ks(isppol,ikpt,ib,ib1)=czero
1149              do ib2=1,oper1%mbandc
1150                oper3%ks(isppol,ikpt,ib,ib1)=oper3%ks(isppol,ikpt,ib,ib1)+ oper1%ks(isppol,ikpt,ib,ib2)*oper2%ks(isppol,ikpt,ib2,ib1)
1151              enddo
1152            enddo
1153          enddo
1154        enddo
1155      enddo
1156    endif
1157  endif
1158 
1159  DBG_EXIT("COLL")
1160 end subroutine prod_oper

m_oper/trace_oper [ Functions ]

[ Top ] [ m_oper ] [ Functions ]

NAME

 trace_oper

FUNCTION

INPUTS

OUTPUT

PARENTS

      impurity_solve,m_green

CHILDREN

      prod_matlu

SOURCE

1036 subroutine trace_oper(oper,trace_ks,trace_loc,opt_ksloc)
1037 
1038  use defs_basis
1039  use m_matlu, only : trace_matlu
1040  use m_errors
1041 
1042 !This section has been created automatically by the script Abilint (TD).
1043 !Do not modify the following lines by hand.
1044 #undef ABI_FUNC
1045 #define ABI_FUNC 'trace_oper'
1046 !End of the abilint section
1047 
1048  implicit none
1049 
1050 !Arguments ------------------------------------
1051 !type
1052  type(oper_type),intent(in) :: oper
1053  real(dp), intent(out) :: trace_ks  !vz_i
1054  real(dp), intent(inout) :: trace_loc(oper%natom,oper%nsppol+1) !vz_i
1055  integer, intent(in) :: opt_ksloc
1056 
1057 !oper variables-------------------------------
1058  integer ::  ib, ikpt, isppol
1059  character(len=500) :: message
1060  real(dp) :: temp1
1061 ! *********************************************************************
1062 
1063  DBG_ENTER("COLL")
1064  if(((opt_ksloc==1.or.opt_ksloc==3).and.(oper%has_opermatlu==0)).or.&
1065 &   ((opt_ksloc==2.or.opt_ksloc==3).and.(oper%has_operks==0))) then
1066    message = " Options in trace_oper are not coherent with definitions of this operator"
1067    MSG_ERROR(message)
1068  endif
1069 
1070  if(opt_ksloc==1.or.opt_ksloc==3) then
1071    trace_ks=zero
1072    temp1=zero
1073    do isppol=1,oper%nsppol
1074      do ikpt=1,oper%nkpt
1075        do ib=1,oper%mbandc
1076          trace_ks=trace_ks+real(oper%ks(isppol,ikpt,ib,ib)*oper%wtk(ikpt))
1077          temp1=temp1+oper%wtk(ikpt)
1078        enddo
1079      enddo
1080    enddo
1081    if(oper%nsppol==1.and.oper%nspinor==1) trace_ks=two*trace_ks
1082 !   write(std_out,*) "temp1",temp1
1083  endif
1084  if(opt_ksloc==2.or.opt_ksloc==3) then
1085    trace_loc=zero
1086    call trace_matlu(oper%matlu,oper%natom,trace_loc)
1087  endif
1088 
1089  DBG_EXIT("COLL")
1090 end subroutine trace_oper

m_oper/upfold_oper [ Functions ]

[ Top ] [ m_oper ] [ Functions ]

NAME

 upfold_oper

FUNCTION

INPUTS

  psichi=<chi|psi> !

OUTPUT

PARENTS

      m_green

CHILDREN

      prod_matlu

SOURCE

742 subroutine upfold_oper(oper,paw_dmft,option,procb,iproc)
743 
744  use defs_basis
745  use m_paw_dmft, only : paw_dmft_type
746  use m_errors
747 
748 !This section has been created automatically by the script Abilint (TD).
749 !Do not modify the following lines by hand.
750 #undef ABI_FUNC
751 #define ABI_FUNC 'upfold_oper'
752 !End of the abilint section
753 
754  implicit none
755 
756 !Arguments ------------------------------------
757 !type
758  integer, intent(in):: option
759  type(oper_type),intent(inout) :: oper
760  type(paw_dmft_type), intent(in) :: paw_dmft
761  integer, optional, intent(in) :: iproc
762  integer, optional, intent(in) :: procb(oper%nkpt)
763 !oper variables-------------------------------
764  integer :: iatom,ib,ib1,ikpt,ispinor,ispinor1,isppol,im1,im
765  integer :: natom,mbandc,ndim,nkpt,nspinor,nsppol,paral
766  integer, allocatable :: procb2(:)
767  character(len=500) :: message
768 ! *********************************************************************
769  ABI_ALLOCATE(procb2,(oper%nkpt))
770  if(present(procb).and.present(iproc)) then
771    paral=1
772    procb2=procb
773 !   write(6,*) "upfold_oper procb",procb
774 !   write(6,*) "iproc",iproc
775 !   write(6,*) size(procb)
776 !   write(6,*) size(procb2)
777 !   write(6,*) procb2(1),procb2(16)
778  else
779    paral=0
780  endif
781 
782 
783 
784  DBG_ENTER("COLL")
785  if((oper%has_opermatlu==0).or.(oper%has_operks==0)) then
786    message = " Operator is not defined to be used in upfold_oper"
787    MSG_ERROR(message)
788  endif
789  if(option<0) then
790  endif
791  nkpt=oper%nkpt
792  natom=paw_dmft%natom
793  nsppol=paw_dmft%nsppol
794  mbandc=paw_dmft%mbandc
795  nspinor=paw_dmft%nspinor
796 
797  do isppol=1,nsppol
798    do ikpt=1,nkpt
799     if ((paral==1.and.(procb2(ikpt)==iproc)).or.(paral==0)) then
800      do ib=1,mbandc
801        do ib1=1,mbandc
802 !               if(ib==1.and.ib1==3) write(std_out,*) "IKPT=",ikpt
803         oper%ks(isppol,ikpt,ib,ib1)=czero
804 
805          do iatom=1,natom
806            if(oper%matlu(iatom)%lpawu.ne.-1) then
807              ndim=2*oper%matlu(iatom)%lpawu+1
808              do im=1,ndim
809                do im1=1,ndim
810                  do ispinor=1,nspinor
811                    do ispinor1=1,nspinor
812 
813 ! psichi(isppol,ikpt,ib,ispinor,iatom,im)=<\chi_{m,R,ispinor)|\Psi(s,k,nu)>
814                      oper%ks(isppol,ikpt,ib,ib1)= oper%ks(isppol,ikpt,ib,ib1) &
815 &                     + ( paw_dmft%psichi(isppol,ikpt,ib1,ispinor1,iatom,im1)        &
816 &                     * oper%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)    &
817 &                     * conjg(paw_dmft%psichi(isppol,ikpt,ib,ispinor,iatom,im)))
818 !               if(ib==1.and.ib1==3) then
819 !                 write(std_out,*) "im,im1",im,im1
820 !                 write(std_out,*) "ispinor,ispinor1",ispinor,ispinor1
821 !                 write(std_out,*) "psichi",paw_dmft%psichi(isppol,ikpt,ib1,ispinor1,iatom,im1)
822 !                 write(std_out,*) "psichi 2",paw_dmft%psichi(isppol,ikpt,ib,ispinor,iatom,im1)
823 !                 write(std_out,*) "oper%matlu", oper%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)
824 !               endif
825 
826                    enddo ! ispinor1
827                  enddo ! ispinor
828                enddo ! im1
829              enddo ! im
830            endif
831          enddo ! iatom
832 
833        enddo ! ib
834      enddo ! ib
835     endif
836    enddo ! ikpt
837  enddo ! isppol
838  ABI_DEALLOCATE(procb2)
839 
840  DBG_EXIT("COLL")
841 end subroutine upfold_oper