TABLE OF CONTENTS


ABINIT/dyson [ Functions ]

[ Top ] [ Functions ]

NAME

 dyson

FUNCTION

 Use the Dyson Equation to compute self-energy from green function 

COPYRIGHT

 Copyright (C) 1999-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 .
 For the initials of contributors, see ~abinit/doc/developers/contributors.txt .

INPUTS

  dtset <type(dataset_type)>=all input variables for this dataset
  istep    =  step of iteration for LDA.
  lda_occup
  mpi_enreg=informations about MPI parallelization
  paw_dmft =  data for self-consistent LDA+DMFT calculations.
  opt_weissgreen= 1: compute weiss from green and self
                = 2: compute green from weiss and self
                = 4: compute green from weiss and self without
                     inversion of weiss  (weiss is previously inverted)

OUTPUT

  edmft  = energy in DMFT.
  paw_dmft =  data for self-consistent LDA+DMFT calculations.

NOTES

PARENTS

      dmft_solve,spectral_function

CHILDREN

      add_matlu,copy_green,destroy_green,init_green,inverse_oper,timab,wrtout

SOURCE

 41 #if defined HAVE_CONFIG_H
 42 #include "config.h"
 43 #endif
 44 
 45 
 46 #include "abi_common.h"
 47 
 48 subroutine dyson(green,paw_dmft,self,weiss,opt_weissself)
 49 
 50  use defs_basis
 51  use m_profiling_abi
 52  use m_errors
 53 
 54  use m_paw_dmft, only: paw_dmft_type
 55  use m_crystal, only : crystal_t
 56  use m_green, only : green_type, destroy_green,init_green,copy_green
 57  use m_oper, only : oper_type,inverse_oper
 58  use m_matlu, only : matlu_type,add_matlu,print_matlu
 59  use m_self, only : self_type
 60 
 61 !This section has been created automatically by the script Abilint (TD).
 62 !Do not modify the following lines by hand.
 63 #undef ABI_FUNC
 64 #define ABI_FUNC 'dyson'
 65  use interfaces_14_hidewrite
 66  use interfaces_18_timing
 67 !End of the abilint section
 68 
 69  implicit none
 70 
 71 !Arguments ------------------------------------
 72 !scalars
 73  type(green_type),intent(inout)  :: green
 74  type(paw_dmft_type), intent(inout) :: paw_dmft
 75  type(self_type),intent(inout)  :: self
 76  type(green_type),intent(inout)  :: weiss
 77  integer,intent(in) :: opt_weissself
 78 ! type(paw_dmft_type), intent(inout)  :: paw_dmft
 79 
 80 !Local variables ------------------------------
 81 !arrays
 82  real(dp) :: tsec(2)
 83 !scalars
 84  integer :: ifreq,natom,nspinor,nsppol,weissinv
 85  type(green_type)  :: greeninv
 86  character(len=500) :: message
 87 ! type
 88 ! type(matlu_type), pointer :: matlutemp,matlu1,matlu2
 89 !************************************************************************
 90 
 91  call timab(623,1,tsec)
 92  DBG_ENTER("COLL")
 93  natom=green%oper(1)%natom
 94  nsppol=green%oper(1)%nsppol
 95  nspinor=green%oper(1)%nspinor
 96  weissinv=1
 97  if(opt_weissself==1) then
 98    write(message,'(2a,i3,13x,a)') ch10,'  ===  Use Dyson Equation => weiss '
 99    call wrtout(std_out,message,'COLL')
100  else if(opt_weissself==2) then
101    write(message,'(2a,i3,13x,a)') ch10,'  ===  Use Dyson Equation => self'
102    call wrtout(std_out,message,'COLL')
103  end if
104 
105 !call xmpi_barrier(spaceComm)
106  if(paw_dmft%dmft_solv==2) weissinv=0
107  call init_green(greeninv,paw_dmft,opt_oper_ksloc=2,wtype=green%w_type)
108  call copy_green(green,greeninv,opt_tw=2)
109 
110  do ifreq=1,green%nw
111    if(opt_weissself==1) then
112      call inverse_oper(greeninv%oper(ifreq),option=1,prtopt=1)
113 !    warning green is now inversed
114      call add_matlu(greeninv%oper(ifreq)%matlu,self%oper(ifreq)%matlu,&
115 &     weiss%oper(ifreq)%matlu,natom,sign_matlu2=1)
116      call inverse_oper(weiss%oper(ifreq),option=1,prtopt=1)
117    else if(opt_weissself==2) then
118 
119    ! write(59,*) paw_dmft%omega_lo(ifreq), real(green%oper(ifreq)%matlu(1)%mat(1,1,1,1,1)),aimag(green%oper(ifreq)%matlu(1)%mat(1,1,1,1,1))
120    ! write(61,*) paw_dmft%omega_lo(ifreq), real(greeninv%oper(ifreq)%matlu(1)%mat(1,1,1,1,1)),aimag(greeninv%oper(ifreq)%matlu(1)%mat(1,1,1,1,1))
121    ! write(60,*) paw_dmft%omega_lo(ifreq), real(weiss%oper(ifreq)%matlu(1)%mat(1,1,1,1,1)),aimag(weiss%oper(ifreq)%matlu(1)%mat(1,1,1,1,1))
122 !    call inverse_oper(weiss%oper(ifreq),option=1,prtopt=1)
123 !    write(62,*) paw_dmft%omega_lo(ifreq), real(weiss%oper(ifreq)%matlu(1)%mat(1,1,1,1,1)),aimag(weiss%oper(ifreq)%matlu(1)%mat(1,1,1,1,1))
124 !    call inverse_oper(weiss%oper(ifreq),option=1,prtopt=1)
125 !    write(63,*) paw_dmft%omega_lo(ifreq), real(self%oper(ifreq)%matlu(1)%mat(1,1,1,1,1)),aimag(self%oper(ifreq)%matlu(1)%mat(1,1,1,1,1))
126 
127 !    write(std_out,*) "-----------------------IFREQ",ifreq
128 !    call print_matlu(greeninv%oper(ifreq)%matlu,paw_dmft%natom,1,opt_diag=-1)
129      call inverse_oper(greeninv%oper(ifreq),option=1,prtopt=1)
130 !    call print_matlu(greeninv%oper(ifreq)%matlu,paw_dmft%natom,1,opt_diag=-1)
131 !    If paw_dmft%dmft_solv==2, then inverse of weiss function is
132 !    computed in hubbard_one.F90
133      if(weissinv/=0) then
134        call inverse_oper(weiss%oper(ifreq),option=1,prtopt=1)
135      end if
136 
137 !    write(std_out,*) weiss%oper(1)%matlu(ifreq)%mat(1,1,1,1,1),"-",greeninv%oper(ifreq)
138      call add_matlu(weiss%oper(ifreq)%matlu,greeninv%oper(ifreq)%matlu,&
139 &     self%oper(ifreq)%matlu,natom,sign_matlu2=-1)
140 
141    ! write(64,*) paw_dmft%omega_lo(ifreq), real(greeninv%oper(ifreq)%matlu(1)%mat(1,1,1,1,1)),aimag(greeninv%oper(ifreq)%matlu(1)%mat(1,1,1,1,1))
142    ! write(65,*) paw_dmft%omega_lo(ifreq), real(weiss%oper(ifreq)%matlu(1)%mat(1,1,1,1,1)),aimag(weiss%oper(ifreq)%matlu(1)%mat(1,1,1,1,1))
143    ! write(66,*) paw_dmft%omega_lo(ifreq), real(self%oper(ifreq)%matlu(1)%mat(1,1,1,1,1)),aimag(self%oper(ifreq)%matlu(1)%mat(1,1,1,1,1))
144    else
145      message = " BUG in dyson.F90"
146      MSG_BUG(message)
147    end if
148  end do
149 
150  call destroy_green(greeninv)
151 
152 
153  call timab(623,2,tsec)
154  DBG_EXIT("COLL")
155 
156 end subroutine dyson