TABLE OF CONTENTS


ABINIT/m_oper [ Modules ]

[ Top ] [ Modules ]

NAME

  m_oper

FUNCTION

COPYRIGHT

 Copyright (C) 2006-2022 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

SOURCE

19 #if defined HAVE_CONFIG_H
20 #include "config.h"
21 #endif
22 
23 
24 #include "abi_common.h"
25 
26 MODULE m_oper
27 
28  use defs_basis
29  use m_abicore
30  use m_errors
31 
32  use m_matlu,    only : matlu_type
33  use m_hide_lapack,  only : xginv
34 
35  implicit none
36 
37  private
38 
39  public :: init_oper
40  public :: diff_oper
41  public :: destroy_oper
42  public :: print_oper
43  public :: inverse_oper
44  public :: loc_oper
45  public :: identity_oper
46  public :: copy_oper
47  public :: trace_oper
48  public :: upfold_oper
49  public :: prod_oper

m_oper/copy_oper [ Functions ]

[ Top ] [ m_oper ] [ Functions ]

NAME

 copy_oper

FUNCTION

INPUTS

OUTPUT

SOURCE

261 subroutine copy_oper(oper1,oper2)
262 
263  use defs_basis
264  use m_matlu, only : copy_matlu
265  use m_errors
266 
267 !Arguments ------------------------------------
268 !type
269  type(oper_type),intent(in) :: oper1
270  type(oper_type),intent(inout) :: oper2 !vz_i
271 
272 !oper variables-------------------------------
273  integer ::  ib, ib1, ikpt, isppol
274 ! *********************************************************************
275  DBG_ENTER("COLL")
276  if(oper1%has_opermatlu==1.and.oper2%has_opermatlu==1)  then
277    call copy_matlu(oper1%matlu,oper2%matlu,oper1%natom)
278  endif
279 
280  if(oper1%has_operks==1.and.oper2%has_operks==1) then
281    do isppol=1,oper1%nsppol
282      do ikpt=1,oper1%nkpt
283        do ib=1,oper1%mbandc
284          do ib1=1,oper1%mbandc
285            oper2%ks(isppol,ikpt,ib,ib1)=oper1%ks(isppol,ikpt,ib,ib1)
286          enddo
287        enddo
288      enddo
289    enddo
290  endif
291 
292  DBG_EXIT("COLL")
293 end subroutine copy_oper

m_oper/destroy_oper [ Functions ]

[ Top ] [ m_oper ] [ Functions ]

NAME

 destroy_oper

FUNCTION

  deallocate oper

INPUTS

  oper

OUTPUT

SOURCE

213 subroutine destroy_oper(oper)
214 
215  use defs_basis
216  use m_crystal, only : crystal_t
217  use m_matlu, only : destroy_matlu
218  use m_errors
219 
220 !Arguments ------------------------------------
221 !scalars
222  type(oper_type),intent(inout) :: oper
223 !local variables-------------------------------
224  character(len=500) :: message
225 
226 !! *********************************************************************
227  DBG_ENTER("COLL")
228  if(oper%has_opermatlu==1) then
229    call destroy_matlu(oper%matlu,oper%natom)
230  else
231    message = " Operator is not defined to be used in destroy_oper"
232    ABI_ERROR(message)
233  endif
234  if ( allocated(oper%matlu))  then
235    ABI_FREE(oper%matlu)
236    oper%has_opermatlu=0
237  endif
238  if ( allocated(oper%ks)) then
239    ABI_FREE(oper%ks)
240    oper%has_operks=0
241  endif
242  oper%wtk => null()
243 !  no deallocation for wtk: wtk is an explicit pointer
244 
245  DBG_EXIT("COLL")
246 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

SOURCE

833 subroutine diff_oper(char1,char2,occup1,occup2,option,toldiff)
834 
835  use defs_basis
836  use m_paw_dmft, only : paw_dmft_type
837  use m_crystal, only : crystal_t
838  use m_matlu, only : diff_matlu
839  use m_errors
840 
841 !Arguments ------------------------------------
842 !type
843  type(oper_type), intent(in) :: occup1,occup2
844  integer, intent(in) :: option
845  real(dp), intent(in) :: toldiff
846  character(len=*), intent(in) :: char1,char2
847 !local variables-------------------------------
848  integer :: mbandc,nkpt
849  character(len=500) :: message
850 ! *********************************************************************
851 
852  DBG_ENTER("COLL")
853  if(occup1%has_opermatlu==0.or.occup2%has_opermatlu==0) then
854    message = " Operators are not defined to be used in diff_oper"
855    ABI_ERROR(message)
856  endif
857  mbandc   = occup1%mbandc
858  nkpt    = occup1%nkpt
859  if(occup1%nkpt/=occup2%nkpt) then
860   write(message,'(a,2x,2i9)')' Operators are not equals',occup1%nkpt,occup2%nkpt
861   ABI_ERROR(message)
862  endif
863 
864  call diff_matlu(char1,char2,occup1%matlu,&
865 & occup2%matlu,occup1%natom,option,toldiff)
866 ! if(option==1) then
867 !  toldiff=tol4
868 !  if( matludiff < toldiff ) then
869 !   write(message,'(6a,e12.4,a,e12.4)') ch10,&
870 !&   '   Differences between ',trim(char1),' and ',trim(char2),' is small enough:',&
871 !&   matludiff,'is lower than',toldiff
872 !   call wrtout(std_out,message,'COLL')
873 !  else
874 !   write(message,'(6a,e12.4,a,e12.4)') ch10,&
875 !&   '   Error: Differences between ',trim(char1),' and ',trim(char2),' is too large:',&
876 !&   matludiff,'is large than',toldiff
877 !   call wrtout(std_out,message,'COLL')
878 !   call abi_abort('COLL')
879 !  endif
880 ! endif
881 ! call abi_abort('COLL')
882 
883  DBG_EXIT("COLL")
884 end subroutine diff_oper

m_oper/identity_oper [ Functions ]

[ Top ] [ m_oper ] [ Functions ]

NAME

 identity_oper

FUNCTION

INPUTS

OUTPUT

SOURCE

753 subroutine identity_oper(oper,option)
754 
755  use defs_basis
756  use m_crystal, only : crystal_t
757  use m_paw_dmft, only : paw_dmft_type
758  use m_errors
759 
760 !Arguments ------------------------------------
761 !type
762  integer, intent(in):: option
763  type(oper_type),intent(inout) :: oper
764 !oper variables-------------------------------
765  integer :: iatom,ib,ikpt,ispinor,isppol,im
766  integer :: natom,mbandc,ndim,nkpt,nspinor,nsppol
767  character(len=500) :: message
768 ! *********************************************************************
769 
770  DBG_ENTER("COLL")
771 
772  if(((option==1.or.option==3).and.(oper%has_opermatlu==0)).or.&
773 &   ((option==2.or.option==3).and.(oper%has_operks==0))) then
774    message = " Options in identity_oper are not coherent with definitions of this operator"
775    ABI_ERROR(message)
776  endif
777  nkpt=oper%nkpt
778  nsppol=oper%nsppol
779  mbandc=oper%mbandc
780  natom=oper%natom
781  nspinor=oper%nspinor
782  oper%ks=czero
783  do iatom=1,natom
784   oper%matlu(iatom)%mat= czero
785  enddo
786 
787  if(option==1.or.option==3) then
788    do isppol=1,nsppol
789      do ikpt=1,nkpt
790        do ib=1,mbandc
791          oper%ks(isppol,ikpt,ib,ib)=cone
792        enddo ! ib
793      enddo ! ikpt
794    enddo ! isppol
795 
796  else if (option==2.or.option==3) then
797    do iatom=1,natom
798      if(oper%matlu(iatom)%lpawu.ne.-1) then
799        ndim=2*oper%matlu(iatom)%lpawu+1
800        do isppol=1,nsppol
801          do im=1,ndim
802            do ispinor=1,nspinor
803              oper%matlu(iatom)%mat(im,im,isppol,ispinor,ispinor)= cone
804            enddo ! ispinor
805          enddo ! im
806        enddo
807      endif ! lpawu
808    enddo ! iatom
809  endif
810 
811  DBG_EXIT("COLL")
812 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

SOURCE

120 subroutine init_oper(paw_dmft,oper,nkpt,wtk,opt_ksloc)
121 
122  use defs_basis
123  use m_matlu, only : init_matlu
124  use m_paw_dmft, only : paw_dmft_type
125  use m_errors
126 
127 !Arguments ------------------------------------
128 !scalars
129  integer, optional, intent(in) :: nkpt
130  integer, optional, intent(in) :: opt_ksloc
131 !type
132  type(paw_dmft_type), intent(in) :: paw_dmft
133  type(oper_type), intent(inout) :: oper
134 !arrays
135  real(dp), pointer, optional :: wtk(:)
136 !oper variables ------------------------------------
137  integer :: iatom,optksloc
138 
139 !************************************************************************
140  DBG_ENTER("COLL")
141  if(present(opt_ksloc)) then
142    optksloc=opt_ksloc
143  else
144    optksloc=3
145  endif
146 
147  if(optksloc/=3) then
148     ! FIXME: empty line!
149  endif
150 
151  oper%has_operks=0
152  oper%has_opermatlu=0
153 ! ===================
154 !  Integers
155 ! ===================
156  oper%nsppol=paw_dmft%nsppol
157  oper%nspinor=paw_dmft%nspinor
158  oper%mbandc=paw_dmft%mbandc
159  oper%natom=paw_dmft%natom
160 
161 ! ===================
162 !  KS variables
163 ! ===================
164  if(optksloc==1.or.optksloc==3) then
165    if(.not.present(nkpt)) then
166      oper%nkpt=paw_dmft%nkpt
167    else
168      oper%nkpt=nkpt
169    endif
170 ! allocate(oper%wtk(oper%nkpt))
171    if(.not.present(wtk)) then
172      oper%wtk=>paw_dmft%wtk
173    else
174      oper%wtk=>wtk
175    endif
176    ABI_MALLOC(oper%ks,(paw_dmft%nsppol,oper%nkpt,paw_dmft%mbandc,paw_dmft%mbandc))
177    oper%has_operks=1
178    oper%ks=czero
179  endif
180 
181 ! ===================
182 !  matlu variables
183 ! ===================
184  if(optksloc==2.or.optksloc==3) then
185    oper%has_opermatlu=0
186    ABI_MALLOC(oper%matlu,(oper%natom))
187    oper%has_opermatlu=1
188    call init_matlu(oper%natom,paw_dmft%nspinor,paw_dmft%nsppol,paw_dmft%lpawu,oper%matlu)
189    do iatom=1,oper%natom
190     oper%matlu(iatom)%mat=czero
191    enddo
192  endif
193 
194  DBG_EXIT("COLL")
195 
196 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

SOURCE

455 subroutine inverse_oper(oper,option,prtopt,procb,iproc)
456 
457  use defs_basis
458  use m_crystal, only : crystal_t
459  use m_paw_dmft, only : paw_dmft_type
460  use m_matlu, only : inverse_matlu
461  use m_errors
462 
463 !Arguments ------------------------------------
464 !type
465  integer, intent(in):: option
466  integer, intent(in) :: prtopt
467  type(oper_type),intent(inout) :: oper
468  integer, optional, intent(in) ::  iproc
469  integer, optional, intent(in) :: procb(oper%nkpt)
470 !oper variables-------------------------------
471  integer :: ikpt,isppol
472  complex(dpc), allocatable :: matrix(:,:)
473  character(len=500) :: message
474  integer :: paral
475  integer, allocatable :: procb2(:)
476 !todo_ba: prb with gwpc here: necessary for matcginv but should be dpc
477 ! *********************************************************************
478  DBG_ENTER("COLL")
479  if(option==2.or.option==3) then
480    ABI_MALLOC(procb2,(oper%nkpt))
481  endif
482 !  if option=1 do inversion in local space
483 !  if option=2 do inversion in KS band space
484 !  if option=3 do both
485  if(present(procb).and.present(iproc)) then
486    paral=1
487    procb2=procb
488  else
489    paral=0
490  endif
491 
492  if(((option==1.or.option==3).and.(oper%has_opermatlu==0)).or.&
493 &   ((option==2.or.option==3).and.(oper%has_operks==0))) then
494    message = " Options are not coherent with definitions of this operator"
495    ABI_ERROR(message)
496  endif
497 
498  if(option==1.or.option==3) then
499    call inverse_matlu(oper%matlu,oper%natom,prtopt)
500  else if(option==2.or.option==3) then
501    ABI_MALLOC(matrix,(oper%mbandc,oper%mbandc))
502      do isppol=1,oper%nsppol
503        do ikpt=1,oper%nkpt
504         if ((paral==1.and.(procb2(ikpt)==iproc)).or.(paral==0)) then
505           matrix(:,:)=oper%ks(isppol,ikpt,:,:)
506 !          write(std_out,*) "isppol,ikpt",isppol,ikpt,m
507 !          write(std_out,*) "isppol,ikpt",matrix
508          !call matcginv_dpc(matrix,oper%mbandc,oper%mbandc)
509          call xginv(matrix,oper%mbandc)
510          oper%ks(isppol,ikpt,:,:)=matrix(:,:)
511         endif
512        enddo ! ikpt
513      enddo ! isppol
514    ABI_FREE(matrix)
515  endif
516 
517  if(option==2.or.option==3) then
518    ABI_FREE(procb2)
519  endif
520  DBG_EXIT("COLL")
521 end subroutine inverse_oper

m_oper/loc_oper [ Functions ]

[ Top ] [ m_oper ] [ Functions ]

NAME

 loc_oper

FUNCTION

INPUTS

OUTPUT

SOURCE

536 subroutine loc_oper(oper,paw_dmft,option,jkpt,procb,iproc)
537 
538  use defs_basis
539  use m_paw_dmft, only : paw_dmft_type
540  use m_errors
541 
542 !Arguments ------------------------------------
543 !type
544  integer, intent(in):: option
545  integer, optional, intent(in):: jkpt
546  type(oper_type),intent(inout) :: oper
547  type(paw_dmft_type), intent(in) :: paw_dmft
548  integer, optional, intent(in) ::  iproc
549  integer, optional, intent(in) :: procb(oper%nkpt)
550 !oper variables-------------------------------
551  integer :: iatom,ib,ib1,ikpt,ikpt1,ispinor,ispinor1,isppol,im1,im
552  integer :: natom,mbandc,ndim,nkpt,nspinor,nsppol,paral
553  character(len=500) :: message
554  integer, allocatable :: procb2(:)
555  logical lvz  !vz_d
556 ! *********************************************************************
557  DBG_ENTER("COLL")
558  ABI_MALLOC(procb2,(oper%nkpt))
559  if((oper%has_opermatlu==0).or.(oper%has_operks==0)) then
560    message = " Operator is not defined to be used in loc_oper"
561    ABI_ERROR(message)
562  endif
563  if(present(procb).and.present(iproc)) then
564    paral=1
565    procb2=procb
566  else
567    paral=0
568  endif
569 
570  if(option<0) then
571  endif
572  nkpt=oper%nkpt
573  natom=oper%natom
574  nsppol=oper%nsppol
575  mbandc=oper%mbandc
576  nspinor=oper%nspinor
577 
578  do iatom=1,natom
579    oper%matlu(iatom)%mat=czero
580  enddo
581 
582  do isppol=1,nsppol
583   do ikpt=1,nkpt
584    ikpt1=ikpt
585    if(present(jkpt)) ikpt1=jkpt
586    lvz=paral==0  !vz_d
587    if(present(iproc)) lvz=lvz.or.(paral==1.and.(procb2(ikpt1)==iproc))  !vz_d
588 !  if ((paral==1.and.(procb2(ikpt1)==iproc)).or.(paral==0)) then    !vz_d
589    if(lvz) then !vz_d
590    do ib=1,mbandc
591     do ib1=1,mbandc
592      do iatom=1,natom
593       if(oper%matlu(iatom)%lpawu.ne.-1) then
594       ndim=2*oper%matlu(iatom)%lpawu+1
595        do im=1,ndim
596         do im1=1,ndim
597          do ispinor=1,nspinor
598           do ispinor1=1,nspinor
599             if (im1 == im .and. im1 == 1 .and. ib1 == ib .and. iatom == 1) then
600 
601             end if
602             oper%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)=     &
603 &            oper%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)+    &
604 &            paw_dmft%psichi(isppol,ikpt1,ib,ispinor,iatom,im)*        &
605 &            conjg(paw_dmft%psichi(isppol,ikpt1,ib1,ispinor1,iatom,im1))* &
606 !false&            paw_dmft%psichi(isppol,ikpt1,ib1,ispinor1,iatom,im1)*        &
607 !false&            conjg(paw_dmft%psichi(isppol,ikpt1,ib,ispinor,iatom,im))* &
608 &            oper%ks(isppol,ikpt,ib,ib1)*oper%wtk(ikpt)
609 ! one  could suppress wtk here if present(jkpt)
610 ! ks(ib,ib1)=ks(ib1,ib) -> ib and ib1 can be underchanged !
611           enddo ! ispinor1
612          enddo ! ispinor
613         enddo ! im1
614        enddo ! im
615       endif
616      enddo ! iatom
617     enddo ! ib
618    enddo ! ib
619    endif
620   enddo ! ikpt
621  enddo ! isppol
622  ABI_FREE(procb2)
623 
624 
625 
626 
627  DBG_EXIT("COLL")
628 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

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

SOURCE

318 subroutine print_oper(oper,option,paw_dmft,prtopt)
319 
320  use defs_basis
321  use m_matlu, only : print_matlu
322  use m_paw_dmft, only : paw_dmft_type
323  use m_errors
324 
325 !Arguments ------------------------------------
326 !type
327  type(paw_dmft_type), intent(in) :: paw_dmft
328  type(oper_type),intent(in) :: oper
329  integer, intent(in) :: option,prtopt
330 
331 !oper variables-------------------------------
332  integer :: ib,ib1,ikpt,isppol,nkptr,iband1,iband2
333  character(len=2000) :: message
334  logical :: ximag
335  real(dp) :: maximag(3)
336 ! *********************************************************************
337  DBG_ENTER("COLL")
338  maximag(:)=zero
339 
340  if(oper%has_opermatlu==1) then
341    write(message,'(2a)') ch10,'   = In the atomic basis'
342    call wrtout(std_out,message,'COLL')
343    call print_matlu(oper%matlu,oper%natom,prtopt)
344  endif
345  if(oper%has_operks==1) then
346    write(message,'(2a)') ch10,'   = In the KS basis'
347    call wrtout(std_out,message,'COLL')
348 
349 !todo_ba complete print_out
350    iband1=1
351    iband2=oper%mbandc
352 !   do ib=1,oper%mbandc
353 !     if(-(paw_dmft%eigen_dft(1,1,ib)+paw_dmft%fermie).ge.0.3) iband1=ib
354 !     if( (paw_dmft%eigen_dft(1,1,ib)-paw_dmft%fermie).le.0.3) iband2=ib
355 !   enddo
356 
357    ximag=.false.
358    if((abs(prtopt)>=3.and.((option<5).or.(option>8))).and.oper%has_operks==1) then
359 !     write(message,'(x,a,a,i4,2x,a)') ch10,'  -KS states'
360 !     call wrtout(std_out,message,'COLL')
361      do isppol=1,paw_dmft%nsppol
362        write(message, '(a,3x,a,i4)') ch10,"--isppol--",isppol
363        call wrtout(std_out,message,'COLL')
364        nkptr=oper%nkpt
365        nkptr=min(oper%nkpt,4)
366        write(message, '(2a)') ch10,&
367 &       "   - ( in the following only the value for the first k-points are printed)"
368        call wrtout(std_out,message,'COLL')
369        do ikpt=1,nkptr
370          if(option<5) then
371            write(message, '(2a,i4,2x,f14.5,a)') ch10,&
372 &           "   -k-pt--",ikpt,oper%wtk(ikpt),"(<-weight(k-pt))"
373            call wrtout(std_out,message,'COLL')
374          else if(abs(prtopt)>=4.or.option>8) then
375            write(message, '(2a,i5,a,i5,a,i5)') ch10,"  Writes occupations for k-pt",&
376 &           ikpt, "and between bands",&
377 &           iband1," and",iband2
378            call wrtout(std_out,message,'COLL')
379          endif
380          do ib=1,oper%mbandc
381            if(option<5) then
382              if(abs(aimag(oper%ks(isppol,ikpt,ib,ib))).ge.tol10) then
383                write(message, '(a,i4,e14.5,3x,e14.5,3x,e21.14)') "   -iband--",ib,&
384 &               paw_dmft%eigen_dft(isppol,ikpt,ib),oper%ks(isppol,ikpt,ib,ib)
385                call wrtout(std_out,message,'COLL')
386              else
387                write(message, '(a,i4,e14.5,3x,e14.5)') "   -iband--",ib,&
388 &               paw_dmft%eigen_dft(isppol,ikpt,ib),real(oper%ks(isppol,ikpt,ib,ib))
389                call wrtout(std_out,message,'COLL')
390              endif
391            endif
392            if(abs(prtopt)>=4.or.option>8.and.ib>=iband1.and.ib<=iband2) then
393              write(message, '(2000(f8.3))') &
394 &               (real(oper%ks(isppol,ikpt,ib,ib1)),ib1=iband1,iband2)
395              call wrtout(std_out,message,'COLL')
396              write(message, '(2000(f8.3))') &
397 &               (aimag(oper%ks(isppol,ikpt,ib,ib1)),ib1=iband1,iband2)
398              call wrtout(std_out,message,'COLL')
399              write(message, '(2000(f8.3))') &
400 &               (sqrt(real(oper%ks(isppol,ikpt,ib,ib1))**2+aimag(oper%ks(isppol,ikpt,ib,ib1))**2),ib1=iband1,iband2)
401              call wrtout(std_out,message,'COLL')
402 !   to write imaginary part
403 !             write(message, '(1000(2f9.3,2x))') &
404 !&               (real(oper%ks(isppol,ikpt,ib,ib1)),imag(oper%ks(isppol,ikpt,ib,ib1)),ib1=iband1,iband2)
405 !             call wrtout(std_out,message,'COLL')
406            endif ! prtopt>=20
407            do ib1=1,oper%mbandc
408              if(abs(aimag(oper%ks(isppol,ikpt,ib,ib1)))>max(tol10,maximag(1))) then
409                ximag=.true.
410                maximag(1)=aimag(oper%ks(isppol,ikpt,ib,ib1))
411              endif
412            enddo
413          enddo ! ib
414        enddo ! ikpt
415      enddo ! isppol
416    else
417     write(message,'(5x,a,i10,a)') '(not written)'
418     call wrtout(std_out,message,'COLL')
419    endif
420    if(ximag) then
421      write(message, '(3a,e12.4,a)')"Occupations are imaginary !",ch10, &
422 &    "  Maximal value is ", maximag(1), ch10
423      ABI_WARNING(message)
424    endif
425  else if(abs(prtopt)>=3.and.((option<5).or.(option>8))) then
426    write(message, '(2a)') ch10," Prb with options and has_operks in print_oper"
427    call wrtout(std_out,message,'COLL')
428  endif ! if oper%has_operks
429 ! write(message, '(2a)') ch10," end print_oper"
430 !     call wrtout(std_out,message,'COLL')
431 
432 
433  DBG_EXIT("COLL")
434 end subroutine print_oper

m_oper/prod_oper [ Functions ]

[ Top ] [ m_oper ] [ Functions ]

NAME

 prod_oper

FUNCTION

INPUTS

OUTPUT

SOURCE

 960 subroutine prod_oper(oper1,oper2,oper3,opt_ksloc)
 961 
 962  use defs_basis
 963  use m_errors
 964  use m_matlu, only : prod_matlu
 965 
 966 !Arguments ------------------------------------
 967 !type
 968  type(oper_type),intent(in) :: oper1
 969  type(oper_type),intent(in) :: oper2
 970  type(oper_type),intent(inout) :: oper3
 971  integer :: opt_ksloc
 972 
 973 !oper variables-------------------------------
 974  integer ::  ib, ib1, ib2, ikpt, isppol
 975 ! *********************************************************************
 976  DBG_ENTER("COLL")
 977  if(opt_ksloc==2) then
 978    if(oper1%has_opermatlu==1.and.oper2%has_opermatlu==1.and.oper3%has_opermatlu==1)  then
 979      call prod_matlu(oper1%matlu,oper2%matlu,oper3%matlu,oper1%natom) !
 980    endif
 981  endif
 982 
 983  if(opt_ksloc==1) then
 984    if(oper1%has_operks==1.and.oper2%has_operks==1.and.oper3%has_operks==1) then
 985      do isppol=1,oper1%nsppol
 986        do ikpt=1,oper1%nkpt
 987          do ib=1,oper1%mbandc
 988            do ib1=1,oper1%mbandc
 989              oper3%ks(isppol,ikpt,ib,ib1)=czero
 990              do ib2=1,oper1%mbandc
 991                oper3%ks(isppol,ikpt,ib,ib1)=oper3%ks(isppol,ikpt,ib,ib1)+ oper1%ks(isppol,ikpt,ib,ib2)*oper2%ks(isppol,ikpt,ib2,ib1)
 992              enddo
 993            enddo
 994          enddo
 995        enddo
 996      enddo
 997    endif
 998  endif
 999 
1000  DBG_EXIT("COLL")
1001 end subroutine prod_oper

m_oper/trace_oper [ Functions ]

[ Top ] [ m_oper ] [ Functions ]

NAME

 trace_oper

FUNCTION

INPUTS

OUTPUT

SOURCE

899 subroutine trace_oper(oper,trace_ks,trace_loc,opt_ksloc)
900 
901  use defs_basis
902  use m_matlu, only : trace_matlu
903  use m_errors
904 
905 !Arguments ------------------------------------
906 !type
907  type(oper_type),intent(in) :: oper
908  real(dp), intent(out) :: trace_ks  !vz_i
909  real(dp), intent(inout) :: trace_loc(oper%natom,oper%nsppol+1) !vz_i
910  integer, intent(in) :: opt_ksloc
911 
912 !oper variables-------------------------------
913  integer ::  ib, ikpt, isppol
914  character(len=500) :: message
915  real(dp) :: temp1
916 ! *********************************************************************
917 
918  DBG_ENTER("COLL")
919  if(((opt_ksloc==1.or.opt_ksloc==3).and.(oper%has_opermatlu==0)).or.&
920 &   ((opt_ksloc==2.or.opt_ksloc==3).and.(oper%has_operks==0))) then
921    message = " Options in trace_oper are not coherent with definitions of this operator"
922    ABI_ERROR(message)
923  endif
924 
925  if(opt_ksloc==1.or.opt_ksloc==3) then
926    trace_ks=zero
927    temp1=zero
928    do isppol=1,oper%nsppol
929      do ikpt=1,oper%nkpt
930        do ib=1,oper%mbandc
931          trace_ks=trace_ks+real(oper%ks(isppol,ikpt,ib,ib)*oper%wtk(ikpt))
932          temp1=temp1+oper%wtk(ikpt)
933        enddo
934      enddo
935    enddo
936    if(oper%nsppol==1.and.oper%nspinor==1) trace_ks=two*trace_ks
937 !   write(std_out,*) "temp1",temp1
938  endif
939  if(opt_ksloc==2.or.opt_ksloc==3) then
940    trace_loc=zero
941    call trace_matlu(oper%matlu,oper%natom,trace_loc)
942  endif
943 
944  DBG_EXIT("COLL")
945 end subroutine trace_oper

m_oper/upfold_oper [ Functions ]

[ Top ] [ m_oper ] [ Functions ]

NAME

 upfold_oper

FUNCTION

INPUTS

  psichi=<chi|psi> !

OUTPUT

SOURCE

644 subroutine upfold_oper(oper,paw_dmft,option,procb,iproc,prt)
645 
646  use defs_basis
647  use m_paw_dmft, only : paw_dmft_type
648  use m_errors
649 
650 !Arguments ------------------------------------
651 !type
652  integer, intent(in):: option
653  type(oper_type),intent(inout) :: oper
654  type(paw_dmft_type), intent(in) :: paw_dmft
655  integer, optional, intent(in) :: iproc
656  integer, optional, intent(in) :: procb(oper%nkpt)
657  integer, optional, intent(in) :: prt
658 !oper variables-------------------------------
659  integer :: iatom,ib,ib1,ikpt,ispinor,ispinor1,isppol,im1,im
660  integer :: natom,mbandc,ndim,nkpt,nspinor,nsppol,paral
661  integer, allocatable :: procb2(:)
662  character(len=500) :: message
663 ! *********************************************************************
664 
665  ABI_UNUSED(prt)
666  ABI_MALLOC(procb2,(oper%nkpt))
667  if(present(procb).and.present(iproc)) then
668    paral=1
669    procb2=procb
670 !   write(6,*) "upfold_oper procb",procb
671 !   write(6,*) "iproc",iproc
672 !   write(6,*) size(procb)
673 !   write(6,*) size(procb2)
674 !   write(6,*) procb2(1),procb2(16)
675  else
676    paral=0
677  endif
678 
679 
680 
681  DBG_ENTER("COLL")
682  if((oper%has_opermatlu==0).or.(oper%has_operks==0)) then
683    message = " Operator is not defined to be used in upfold_oper"
684    ABI_ERROR(message)
685  endif
686  if(option<0) then
687  endif
688  nkpt=oper%nkpt
689  natom=paw_dmft%natom
690  nsppol=paw_dmft%nsppol
691  mbandc=paw_dmft%mbandc
692  nspinor=paw_dmft%nspinor
693 
694  do isppol=1,nsppol
695    do ikpt=1,nkpt
696     if ((paral==1.and.(procb2(ikpt)==iproc)).or.(paral==0)) then
697      do ib=1,mbandc
698        do ib1=1,mbandc
699 !               if(ib==1.and.ib1==3) write(std_out,*) "IKPT=",ikpt
700         oper%ks(isppol,ikpt,ib,ib1)=czero
701 
702          do iatom=1,natom
703            if(oper%matlu(iatom)%lpawu.ne.-1) then
704              ndim=2*oper%matlu(iatom)%lpawu+1
705              do im=1,ndim
706                do im1=1,ndim
707                  do ispinor=1,nspinor
708                    do ispinor1=1,nspinor
709 
710 ! psichi(isppol,ikpt,ib,ispinor,iatom,im)=<\chi_{m,R,ispinor)|\Psi(s,k,nu)>
711                      oper%ks(isppol,ikpt,ib,ib1)= oper%ks(isppol,ikpt,ib,ib1) &
712 &                     + ( paw_dmft%psichi(isppol,ikpt,ib1,ispinor1,iatom,im1)        &
713 &                     * oper%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)    &
714 &                     * conjg(paw_dmft%psichi(isppol,ikpt,ib,ispinor,iatom,im)))
715               ! if(present(prt).and.(ib==1.and.ib1==1)) then
716               !   write(6,*) "im,im1",im,im1
717               !   write(6,*) "ispinor,ispinor1",ispinor,ispinor1
718               !   write(6,*) "psichi",paw_dmft%psichi(isppol,ikpt,ib1,ispinor1,iatom,im1)
719               !   write(6,*) "psichi 2",paw_dmft%psichi(isppol,ikpt,ib,ispinor,iatom,im1)
720               !   write(6,*) "oper%matlu", oper%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)
721               ! endif
722 
723                    enddo ! ispinor1
724                  enddo ! ispinor
725                enddo ! im1
726              enddo ! im
727            endif
728          enddo ! iatom
729 
730        enddo ! ib
731      enddo ! ib
732     endif
733    enddo ! ikpt
734  enddo ! isppol
735  ABI_FREE(procb2)
736 
737  DBG_EXIT("COLL")
738 end subroutine upfold_oper