TABLE OF CONTENTS


ABINIT/m_self [ Modules ]

[ Top ] [ Modules ]

NAME

  m_self

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 .

PARENTS

CHILDREN

SOURCE

20 #if defined HAVE_CONFIG_H
21 #include "config.h"
22 #endif
23 
24 
25 #include "abi_common.h"
26 
27 MODULE m_self
28 
29  use defs_basis
30  use m_xmpi
31  use m_errors
32  use m_abicore
33 
34  use m_oper,     only : oper_type
35  use m_matlu,    only : matlu_type
36  use m_fstrings, only : int2char4
37 
38  implicit none
39 
40  private 
41 
42  public :: alloc_self
43  public :: initialize_self
44  public :: destroy_self
45  public :: print_self
46  public :: rw_self
47  public :: dc_self
48  public :: new_self
49  public :: make_qmcshift_self

m_self/alloc_self [ Functions ]

[ Top ] [ m_self ] [ Functions ]

NAME

 alloc_self

FUNCTION

  Allocate variables used in type self_type.

INPUTS

  self <type(self_type)>= variables related to self-energy
  paw_dmft <type(paw_dmft_type)> =  variables related to self-consistent LDA+DMFT calculations.
  opt_oper = 1  Allocate only quantities in the KS basis.
             2  Allocate only quantities in the local basis.
             3  Allocate quantities in both the KS and local basis.
  wtype = "real" Self energy will be computed for real frequencies
        = "imag" Self energy will be computed for imaginary frequencies

 OUTPUTS
  self <type(self_type)>= variables related to self-energy

PARENTS

      m_self

CHILDREN

      shift_matlu,wrtout

SOURCE

130 subroutine alloc_self(self,paw_dmft,opt_oper,wtype)
131 
132  use defs_basis
133  use defs_abitypes
134  use m_crystal, only : crystal_t
135  use m_oper, only : init_oper
136  use m_paw_dmft, only: paw_dmft_type
137 
138 !This section has been created automatically by the script Abilint (TD).
139 !Do not modify the following lines by hand.
140 #undef ABI_FUNC
141 #define ABI_FUNC 'alloc_self'
142 !End of the abilint section
143 
144  implicit none
145 
146 !Arguments ------------------------------------
147 !scalars
148 !type
149  type(self_type), intent(inout) :: self
150  type(paw_dmft_type), intent(in) :: paw_dmft
151  integer, optional, intent(in) :: opt_oper
152  character(len=4), optional :: wtype
153 !local variables ------------------------------------
154  integer :: ifreq,optoper
155 
156 !************************************************************************
157  if(present(opt_oper)) then
158    optoper=opt_oper
159  else
160    optoper=2
161  endif
162  if(present(wtype)) then
163    self%w_type=wtype
164  else
165    self%w_type="imag"
166  endif
167  if(self%w_type=="imag") then
168    self%nw=paw_dmft%dmft_nwlo
169    self%omega=>paw_dmft%omega_lo
170  else if(self%w_type=="real") then
171    self%nw=2*paw_dmft%dmft_nwr
172    self%omega=>paw_dmft%omega_r
173  endif
174 
175  self%dmft_nwlo=paw_dmft%dmft_nwlo
176  self%dmft_nwli=paw_dmft%dmft_nwli
177  self%iself_cv=0
178  
179  call init_oper(paw_dmft,self%hdc,opt_ksloc=optoper)
180  ABI_DATATYPE_ALLOCATE(self%oper,(self%nw))
181  do ifreq=1,self%nw
182   call init_oper(paw_dmft,self%oper(ifreq),opt_ksloc=optoper)
183  enddo
184 
185  if(paw_dmft%dmft_solv==4) then
186    ABI_ALLOCATE(self%qmc_shift,(paw_dmft%natom))
187    ABI_ALLOCATE(self%qmc_xmu,(paw_dmft%natom))
188    self%qmc_shift(:)=zero
189    self%qmc_xmu(:)=zero
190  endif
191 
192 end subroutine alloc_self

m_self/dc_self [ Functions ]

[ Top ] [ m_self ] [ Functions ]

NAME

 dc_self

FUNCTION

INPUTS

  charge_loc(cryst_struc%natom,paw_dmft%nsppol+1)= total charge for correlated electrons on a given atom, and for spin
  cryst_struc <type(crystal_t)>=variables related to crystal structure 
  hu <type(hu_type)>= variables related to the interaction between electrons
  self <type(self_type)>= variables related to self-energy
  dmft_dc = 1 Full localized Limit double counting.
           2 Around Mean Field (without SO)
           0 not double counting 
  prtopt = integer which precises the amount of printing (not used here)

OUTPUT

  self <type(self_type)>= variables related to self-energy
  hu <type(hu_type)>= variables related to the interaction between electrons

PARENTS

      m_dmft,spectral_function

CHILDREN

      shift_matlu,wrtout

SOURCE

433 subroutine dc_self(charge_loc,cryst_struc,hu,self,dmft_dc,prtopt)
434 
435  use defs_basis
436  use m_crystal, only : crystal_t
437  use m_paw_dmft, only : paw_dmft_type
438  use m_hu, only : hu_type
439 
440 !This section has been created automatically by the script Abilint (TD).
441 !Do not modify the following lines by hand.
442 #undef ABI_FUNC
443 #define ABI_FUNC 'dc_self'
444 !End of the abilint section
445 
446  implicit none
447 
448 !Arguments ------------------------------------
449 !type
450  type(crystal_t),intent(in) :: cryst_struc
451  type(self_type),intent(inout) :: self
452  real(dp), intent(in) :: charge_loc(cryst_struc%natom,self%hdc%nsppol+1)
453  type(hu_type),intent(inout) :: hu(cryst_struc%ntypat)
454  integer, intent(in) :: prtopt,dmft_dc
455 
456 !Local variables-------------------------------
457  integer :: iatom,isppol,ispinor,lpawu,m1,nspinor,nsppol
458  real(dp) :: ntot,jpawu_dc
459  character(len=500) :: message
460 ! *********************************************************************
461  nspinor = self%hdc%nspinor
462  nsppol  = self%hdc%nsppol
463 
464  do iatom=1,self%hdc%natom
465    lpawu=self%hdc%matlu(iatom)%lpawu
466    ntot=charge_loc(iatom,nsppol+1)
467 !     nup=charge_loc(iatom,1)
468 !     ndn=charge_loc(iatom,nsppol)
469    if(lpawu/=-1) then
470      self%hdc%matlu(iatom)%mat=czero
471      do isppol = 1 , nsppol
472        do ispinor = 1 , nspinor
473          do m1=1, 2*lpawu+1
474            if(dmft_dc==1.or.dmft_dc==4.or.dmft_dc==5) then ! FLL
475              if(dmft_dc==1.or.dmft_dc==5) then
476                jpawu_dc=hu(cryst_struc%typat(iatom))%jpawu
477              else if (dmft_dc==4) then ! FLL
478                jpawu_dc=zero
479              endif
480              if(nspinor==2.or.dmft_dc==5) then
481                self%hdc%matlu(iatom)%mat(m1,m1,isppol,ispinor,ispinor) =  &
482 &                hu(cryst_struc%typat(iatom))%upawu * ( ntot - one / two )     &
483 &                - one/two * jpawu_dc * ( ntot - one       )
484              else
485                self%hdc%matlu(iatom)%mat(m1,m1,isppol,ispinor,ispinor)=  &
486          &       hu(cryst_struc%typat(iatom))%upawu * ( ntot - one / two ) &
487 !&     -  hu(cryst_struc%typat(iatom))%jpawu * ( charge_loc(iatom,2-nsppol+1) - one)
488 &                 -  jpawu_dc * ( charge_loc(iatom,isppol) - one / two )
489              endif
490            else if(dmft_dc==2) then  ! AMF
491              if(nspinor==2) then
492                !write(message,'(a,i4,i4,2x,e20.10)') " AMF Double counting not implemented for SO"
493                !MSG_ERROR(message)
494                self%hdc%matlu(iatom)%mat(m1,m1,isppol,ispinor,ispinor)= &
495                hu(cryst_struc%typat(iatom))%upawu * ntot/two &
496                + ( hu(cryst_struc%typat(iatom))%upawu - hu(cryst_struc%typat(iatom))%jpawu )&
497                *ntot/two*(float(2*lpawu))/(float(2*lpawu+1))
498                write(message,'(a,i4,i4,2x,e20.10)') " AMF Double counting is under test for SOC"
499                MSG_WARNING(message)
500              else
501                if(nsppol==2) then
502                  self%hdc%matlu(iatom)%mat(m1,m1,isppol,ispinor,ispinor)=  &
503 &                  hu(cryst_struc%typat(iatom))%upawu * charge_loc(iatom,2-isppol+1) &
504 &                  +  (hu(cryst_struc%typat(iatom))%upawu - hu(cryst_struc%typat(iatom))%jpawu )&
505 &                  *charge_loc(iatom,isppol)*(float(2*lpawu))/(float(2*lpawu+1)) 
506                 else  if(nsppol==1) then
507                   self%hdc%matlu(iatom)%mat(m1,m1,isppol,ispinor,ispinor)=  &
508 &                   hu(cryst_struc%typat(iatom))%upawu * charge_loc(iatom,isppol) &
509 &                    +  (hu(cryst_struc%typat(iatom))%upawu - hu(cryst_struc%typat(iatom))%jpawu )&
510 &                   *charge_loc(iatom,isppol)*(float(2*lpawu))/(float(2*lpawu+1)) 
511                 endif
512 !                 write(std_out,*) "AMF",  charge_loc(iatom,2-isppol+1)
513 !                 write(std_out,*) "AMF",  charge_loc(iatom,isppol+1)
514 !                 write(std_out,*) "AMF",  lpawu
515 !                 write(std_out,*) "AMF",  hu(cryst_struc%typat(iatom))%upawu 
516 !                 write(std_out,*) "AMF", self%hdc%matlu(iatom)%mat(m1,m1,isppol,ispinor,ispinor)
517              endif
518            else
519              MSG_ERROR("not implemented")
520            endif
521          enddo  ! m1
522        enddo  ! ispinor
523      enddo ! isppol
524    endif
525  enddo ! iatom
526 
527  if(prtopt>0) then
528  endif
529 
530 
531 end subroutine dc_self

m_self/destroy_self [ Functions ]

[ Top ] [ m_self ] [ Functions ]

NAME

 destroy_self

FUNCTION

  deallocate self

INPUTS

  self <type(self_type)>= variables related to self-energy

OUTPUT

PARENTS

      m_dmft,spectral_function

CHILDREN

      shift_matlu,wrtout

SOURCE

294 subroutine destroy_self(self)
295 
296  use defs_basis
297  use m_crystal, only : crystal_t
298  use m_oper, only : destroy_oper
299 
300 !This section has been created automatically by the script Abilint (TD).
301 !Do not modify the following lines by hand.
302 #undef ABI_FUNC
303 #define ABI_FUNC 'destroy_self'
304 !End of the abilint section
305 
306  implicit none
307 
308 !Arguments ------------------------------------
309 !scalars
310  type(self_type),intent(inout) :: self
311 
312 !local variables-------------------------------
313  integer :: ifreq
314 
315 ! *********************************************************************
316  if ( allocated(self%oper))       then
317    do ifreq=1,self%nw
318     call destroy_oper(self%oper(ifreq))
319    enddo
320    ABI_DATATYPE_DEALLOCATE(self%oper)
321  end if
322 
323  call destroy_oper(self%hdc)
324  if (allocated(self%qmc_shift)) then
325    ABI_DEALLOCATE(self%qmc_shift)
326  end if
327  if (allocated(self%qmc_xmu))  then
328    ABI_DEALLOCATE(self%qmc_xmu)
329  end if
330  self%omega => null()
331 
332 end subroutine destroy_self

m_self/initialize_self [ Functions ]

[ Top ] [ m_self ] [ Functions ]

NAME

 initialize_self

FUNCTION

  Initialize self-energy.

INPUTS

  cryst_struc <type(crystal_t)>=variables related to crystal structure 
  self <type(self_type)>= variables related to self-energy
  paw_dmft <type(paw_dmft_type)> =  variables related to self-consistent LDA+DMFT calculations.
  opt_read =  not used for the moment 
  wtype = "real" Self energy will be computed for real frequencies
        = "imag" Self energy will be computed for imaginary frequencies

 OUTPUTS
  self <type(self_type)>= variables related to self-energy

PARENTS

      m_dmft,spectral_function

CHILDREN

      shift_matlu,wrtout

SOURCE

222 subroutine initialize_self(self,paw_dmft,wtype)
223 
224  use defs_basis
225  use defs_abitypes
226  use m_crystal, only : crystal_t
227  use m_oper, only : init_oper,loc_oper
228  use m_matlu, only : print_matlu
229  use m_paw_dmft, only: paw_dmft_type
230 
231 !This section has been created automatically by the script Abilint (TD).
232 !Do not modify the following lines by hand.
233 #undef ABI_FUNC
234 #define ABI_FUNC 'initialize_self'
235 !End of the abilint section
236 
237  implicit none
238 
239 !Arguments ------------------------------------
240 !scalars
241 !type
242  type(self_type), intent(inout) :: self
243  type(paw_dmft_type), intent(inout) :: paw_dmft
244  character(len=4), optional, intent(in) :: wtype
245 !local variables ------------------------------------
246 ! character(len=500) :: message
247  integer :: iatom,ifreq
248  character(len=4) :: wtype2
249 !************************************************************************
250  if(present(wtype)) then
251    wtype2=wtype
252  else
253    wtype2="imag"
254  endif
255 
256  
257  call alloc_self(self,paw_dmft,opt_oper=2,wtype=wtype2) !  opt_oper=1 is not useful and not implemented
258  do ifreq=1,self%nw
259    do iatom=1,paw_dmft%natom
260      self%oper(ifreq)%matlu(iatom)%mat=czero
261    enddo
262  enddo
263 ! if(paw_dmft%dmft_rslf==1.and.opt_read==1) then
264 !   call rw_self(cryst_struc,self,mpi_enreg,paw_dmft,pawtab,pawprtvol=2,opt_rw=1)
265 ! endif
266 ! write(message,'(a,a)') ch10,"   Self-energy for large frequency is"
267 ! call wrtout(std_out,message,'COLL')
268 ! call print_matlu(self%oper(paw_dmft%dmft_nwlo)%matlu,  &
269 !&                 paw_dmft%natom,3)
270 
271 end subroutine initialize_self

m_self/make_qmcshift_self [ Functions ]

[ Top ] [ m_self ] [ Functions ]

NAME

 make_qmcshift_hu

FUNCTION

INPUTS

  hu <type(hu_type)> = U interaction
  paw_dmft  <type(paw_dmft_type)> = paw+dmft related data

OUTPUT

  self%qmc_shift in self <type(self_type)> = Self-energy

PARENTS

CHILDREN

      wrtout

SOURCE

1207 subroutine make_qmcshift_self(cryst_struc,hu,self,apply)
1208 
1209  use defs_basis
1210  use m_paw_dmft, only : paw_dmft_type
1211  use m_crystal, only : crystal_t
1212  use m_hu, only : hu_type
1213  use m_matlu, only : shift_matlu
1214 
1215 !This section has been created automatically by the script Abilint (TD).
1216 !Do not modify the following lines by hand.
1217 #undef ABI_FUNC
1218 #define ABI_FUNC 'make_qmcshift_self'
1219 !End of the abilint section
1220 
1221  implicit none
1222 
1223 !Arguments ------------------------------------
1224 !type
1225  type(crystal_t),intent(in) :: cryst_struc
1226  type(hu_type),intent(in) :: hu(cryst_struc%ntypat)
1227  type(self_type),intent(inout) :: self
1228  logical, optional :: apply
1229 
1230 !Local variables-------------------------------
1231  integer :: im,iatom,ifreq,itypat,lpawu,tndim
1232  real(dp) :: hu_shift2
1233  character(len=500) :: message
1234 ! *********************************************************************
1235 
1236  do iatom = 1 , cryst_struc%natom
1237    lpawu=self%hdc%matlu(iatom)%lpawu
1238    tndim=2*lpawu+1
1239    itypat=cryst_struc%typat(iatom)
1240    if(lpawu/=-1) then
1241      self%qmc_shift(iatom) = zero
1242      do im =1, 2*tndim-1
1243 !       write(std_out,*)"make before",self%qmc_shift(iatom)
1244        self%qmc_shift(iatom) = self%qmc_shift(iatom) + hu(itypat)%uqmc(im)
1245 !       write(std_out,*)"make after",self%qmc_shift(iatom)
1246      enddo
1247      self%qmc_shift(iatom) = self%qmc_shift(iatom) / two
1248      hu_shift2 = hu(itypat)%uqmc(1)
1249 
1250      do im = 2*tndim, 2*tndim + 2*tndim -3
1251        hu_shift2 = hu_shift2 + hu(itypat)%uqmc(im)
1252      enddo
1253 
1254      hu_shift2 = hu_shift2 / two
1255      write(message,'(2a,i4)')  ch10,'  -------> For Correlated atom',iatom
1256      call wrtout(std_out,  message,'COLL')
1257 
1258      if(abs(self%qmc_shift(iatom)-hu_shift2)>tol6) then
1259        write(message,'(2a,2f16.7)')  "  Shift for QMC is not correctly"&
1260 &      ," computed",self%qmc_shift(iatom),hu_shift2
1261        MSG_ERROR(message)
1262      endif ! shifts not equals
1263 
1264      write(message,'(4x,a,f16.7)')  &
1265 &     "  Shift for QMC (used to compute G(w)) is (in Ha) :",&
1266 &     self%qmc_shift(iatom)
1267      call wrtout(std_out,  message,'COLL')
1268 
1269      self%qmc_xmu(iatom)=-self%qmc_shift(iatom)
1270      self%qmc_xmu(iatom)=zero
1271      write(message,'(4x,a,f16.7)')  &
1272 &     "Artificial Shift used in QMC AND to compute G is (in Ha) :",self%qmc_xmu(iatom)
1273      MSG_WARNING(message)
1274 
1275    endif ! lpawu/=1
1276  enddo ! natom
1277 
1278  if(present(apply)) then
1279    if (apply) then
1280      write(message,'(5x,a,f16.7,a)')  " Shifts applied to self"
1281      call wrtout(std_out,  message,'COLL')
1282      do ifreq=1,self%nw
1283        call shift_matlu(self%oper(ifreq)%matlu,cryst_struc%natom,cmplx(self%qmc_shift,0.d0,kind=dp),1)
1284        call shift_matlu(self%oper(ifreq)%matlu,cryst_struc%natom,cmplx(self%qmc_xmu,0.d0,kind=dp),1)
1285      enddo
1286    endif
1287  endif
1288 
1289 
1290 end subroutine make_qmcshift_self

m_self/new_self [ Functions ]

[ Top ] [ m_self ] [ Functions ]

NAME

 new_self

FUNCTION

  Mix Old and New self_energy with the mixing coefficient dmft_mxsf

INPUTS

  self <type(self_type)>= variables related to self-energy
  self_new <type(self_type)>= variables related to the new self-energy
  paw_dmft  <type(paw_dmft_type)>= paw+dmft related data
  opt_mix not used

OUTPUT

  self <type(self_type)>= variables related to mixed self-energy

PARENTS

      m_dmft

CHILDREN

      shift_matlu,wrtout

SOURCE

1095 subroutine new_self(self,self_new,paw_dmft,opt_mix)
1096 
1097  use defs_basis
1098  use defs_abitypes
1099  use m_crystal, only : crystal_t
1100  use m_paw_dmft, only : paw_dmft_type
1101  use m_matlu, only : copy_matlu
1102 
1103 !This section has been created automatically by the script Abilint (TD).
1104 !Do not modify the following lines by hand.
1105 #undef ABI_FUNC
1106 #define ABI_FUNC 'new_self'
1107 !End of the abilint section
1108 
1109  implicit none
1110 
1111 !Arguments ------------------------------------
1112 !type
1113 ! type(crystal_t),intent(in) :: cryst_struc
1114  type(self_type),intent(inout) :: self
1115  type(self_type),intent(in) :: self_new
1116  type(paw_dmft_type), intent(in) :: paw_dmft
1117  integer,intent(in) :: opt_mix
1118 
1119 !local variables-------------------------------
1120  integer :: iatom,icount,ifreq,im,im1,ispinor,isppol,ispinor1
1121  integer :: natom,ndim,nspinor,nsppol
1122  real(dp) :: alpha,diff_self,sum_self
1123  complex(dpc) :: s1,s2
1124  character(len=500) :: message
1125 ! *********************************************************************
1126  natom=self%hdc%natom
1127  nsppol=paw_dmft%nsppol
1128  nspinor=paw_dmft%nspinor
1129  alpha=paw_dmft%dmft_mxsf
1130 
1131  if(opt_mix==1) then
1132  endif
1133  sum_self=zero
1134  diff_self=zero
1135  icount=0
1136  do iatom=1,natom
1137    if(self%oper(1)%matlu(iatom)%lpawu.ne.-1) then
1138      ndim=2*self%oper(1)%matlu(iatom)%lpawu+1
1139      do isppol=1,nsppol
1140        do ispinor=1,nspinor
1141          do ispinor1=1,nspinor
1142            do im=1,ndim
1143              do im1=1,ndim
1144                do ifreq=1,self%nw
1145 !  warning: self_new is the recent self-energy, which is mixed with self
1146 !  to give self= mixed self energy. self_new is deallocated just after.
1147                  self%oper(ifreq)%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)=     &
1148 &                  (one-alpha)*self%oper(ifreq)%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1) +    &
1149 &                  (alpha)*self_new%oper(ifreq)%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)
1150                      s1=self%hdc%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)
1151                  s2=self_new%hdc%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)
1152                  if((ispinor==ispinor1).and.(im==im1)) then
1153                    diff_self=diff_self+dsqrt(real(s1-s2)**2+aimag(s1-s2)**2)
1154                    sum_self=sum_self+dsqrt(real(s1)**2+aimag(s1)**2)
1155                    icount=icount+1
1156                  endif
1157                enddo
1158                self%hdc%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)=             &
1159 &               (one-alpha)*self%hdc%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)   +          &
1160 &               (alpha)*self_new%hdc%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)
1161              enddo
1162            enddo
1163          enddo
1164        enddo
1165      enddo ! isppol
1166    endif ! lpawu=/-1
1167  enddo ! iatom
1168  diff_self=diff_self/icount
1169 
1170  write(message,'(8x,a,e12.5)') "DMFT Loop: Precision on self-energy is",diff_self
1171  call wrtout(std_out,message,'COLL')
1172  if(diff_self<paw_dmft%dmft_fepr.and.sum_self>tol6.and.paw_dmft%idmftloop>=2) then
1173     write(message,'(a,8x,a,e9.2,a,8x,a)') ch10, "Change of self =<", paw_dmft%dmft_fepr,&
1174 &    ch10,"DMFT Loop: Self Energy is converged" 
1175     call wrtout(std_out,message,'COLL')
1176     self%iself_cv=1
1177  else
1178     write(message,'(a,8x,a)') ch10,"DMFT Loop: Self Energy is not converged" 
1179     call wrtout(std_out,message,'COLL')
1180     self%iself_cv=0
1181  endif
1182 
1183 
1184 end subroutine new_self

m_self/print_self [ Functions ]

[ Top ] [ m_self ] [ Functions ]

NAME

 print_self

FUNCTION

  print occupations

INPUTS

  self <type(self_type)>= variables related to self-energy
  option = 1 Do not print double counting.
           2 Print double counting
  paw_dmft <type(paw_dmft_type)> =  variables related to self-consistent LDA+DMFT calculations.
  prtopt = integer which precises the amount of printing in the subroutine called 

OUTPUT

  self <type(self_type)>= variables related to self-energy

PARENTS

      m_dmft

CHILDREN

      shift_matlu,wrtout

SOURCE

360 subroutine print_self(self,prtdc,paw_dmft,prtopt)
361 
362  use defs_basis
363  use m_oper, only : print_oper
364  use m_paw_dmft, only : paw_dmft_type
365  use m_matlu, only : print_matlu
366 
367 !This section has been created automatically by the script Abilint (TD).
368 !Do not modify the following lines by hand.
369 #undef ABI_FUNC
370 #define ABI_FUNC 'print_self'
371 !End of the abilint section
372 
373  implicit none
374 
375 !Arguments ------------------------------------
376 !type
377  type(paw_dmft_type), intent(in) :: paw_dmft
378  type(self_type),intent(in) :: self
379  character(len=*), intent(in) :: prtdc
380  integer,intent(in) :: prtopt
381 
382 
383 !local variables-------------------------------
384  character(len=500) :: message
385 ! *********************************************************************
386 
387  write(message,'(2a)') ch10,"  == The self-energy for smallest frequency is   == "
388  call wrtout(std_out,message,'COLL')
389  call print_oper(self%oper(1),1,paw_dmft,prtopt)
390 ! write(message,'(2a)') ch10,"  == The self-energy for small (3) frequency is   == "
391 ! call wrtout(std_out,message,'COLL')
392 ! call print_oper(self%oper(3),1,paw_dmft,prtopt)
393  write(message,'(2a)') ch10,"  == The self-energy for large frequency is   == "
394  call wrtout(std_out,message,'COLL')
395  call print_oper(self%oper(self%nw),1,paw_dmft,prtopt)
396  if(prtdc=="print_dc") then
397    write(message,'(2a)') ch10,"  == The double counting hamiltonian is (diagonal part)  == "
398    call wrtout(std_out,message,'COLL')
399    call print_matlu(self%hdc%matlu,paw_dmft%natom,prtopt,opt_diag=1)
400  endif
401 
402 end subroutine print_self

m_self/rw_self [ Functions ]

[ Top ] [ m_self ] [ Functions ]

NAME

 rw_self

FUNCTION

INPUTS

  self <type(self_type)>= variables related to self-energy
  paw_dmft  <type(paw_dmft_type)>= paw+dmft related data
  prtopt = integer which precises the amount of printing 
  opt_rw = 1  Read Self-Energy.
           2  Write Self-Energy.
           3  Impose Self-Energy.

OUTPUT

PARENTS

      m_dmft,spectral_function

CHILDREN

      shift_matlu,wrtout

SOURCE

 559 subroutine rw_self(self,paw_dmft,prtopt,opt_rw,istep_iter,opt_char)
 560 
 561  use defs_basis
 562  use defs_abitypes
 563 
 564  use m_io_tools, only : get_unit
 565  use m_crystal, only : crystal_t
 566  use m_paw_dmft, only : paw_dmft_type
 567  use m_matlu, only : copy_matlu,shift_matlu
 568 
 569 !This section has been created automatically by the script Abilint (TD).
 570 !Do not modify the following lines by hand.
 571 #undef ABI_FUNC
 572 #define ABI_FUNC 'rw_self'
 573 !End of the abilint section
 574 
 575  implicit none
 576 
 577 !Arguments ------------------------------------
 578 !type
 579  type(self_type),intent(inout) :: self
 580  !type(MPI_type), intent(in) :: mpi_enreg
 581  type(paw_dmft_type), intent(inout) :: paw_dmft
 582  integer,intent(in) :: prtopt
 583  integer,intent(in),optional :: opt_rw,istep_iter
 584  character(len=4), optional :: opt_char
 585 
 586 !local variables-------------------------------
 587  logical :: lexist
 588  complex(dpc), allocatable :: buffer(:)
 589  integer :: iall,iatom,iatu,ier,iexist2,ifreq,im,im1,ioerr,ispinor,ispinor1,isppol,istepiter,istep,istep_imp
 590  integer :: icount,iexit,iter,iter_imp,master,mbandc,myproc,natom,ncount,ndim,nkpt,nproc,nrecl,nspinor,nsppol,spacecomm
 591  integer :: natom_read,nsppol_read,nspinor_read,ndim_read,nw_read,optrw
 592  character(len=30000) :: message ! Big buffer to avoid buffer overflow.
 593  integer,allocatable :: unitselffunc_arr(:)
 594  character(len=fnlen) :: tmpfil
 595  character(len=1) :: tag_is
 596  character(len=10) :: tag_at
 597  character(len=4) :: chtemp
 598  real(dp):: xtemp,fermie_read
 599  real(dp), allocatable:: s_r(:,:,:,:),s_i(:,:,:,:),fermie_read2(:)
 600 ! *********************************************************************
 601 
 602 ! Initialise spaceComm, myproc, and nproc
 603  istep=0 ; iter=0 ; istep_imp=0 ; iter_imp=0
 604  if(present(opt_rw)) then
 605    optrw=opt_rw
 606  else
 607    optrw=0
 608  endif
 609  if(present(istep_iter)) then
 610    istepiter=istep_iter
 611  else
 612    istepiter=0
 613  endif
 614  if(paw_dmft%use_fixed_self>0) then
 615    istep=istepiter/1000
 616    iter=istepiter-(istepiter/1000)*1000
 617    istep_imp=paw_dmft%use_fixed_self/1000
 618    iter_imp=paw_dmft%use_fixed_self-(paw_dmft%use_fixed_self/1000)*1000
 619  endif
 620 
 621  if(paw_dmft%dmft_rslf<=0.and.optrw==1) optrw=0
 622  iexit=0
 623  ioerr=0
 624  iexist2=1
 625  lexist=.true.
 626  spacecomm=paw_dmft%spacecomm
 627  myproc=paw_dmft%myproc
 628  nproc=paw_dmft%nproc
 629  master=0
 630 
 631 ! write(std_out,*) "myproc,master",myproc,master
 632  if(prtopt>200) then
 633  endif
 634  natom=self%oper(1)%natom
 635  nsppol=paw_dmft%nsppol
 636  nspinor=paw_dmft%nspinor
 637  mbandc=paw_dmft%mbandc
 638  nkpt=paw_dmft%nkpt
 639  if((optrw==2.or.optrw==1).and.myproc==master)then
 640    ABI_ALLOCATE(unitselffunc_arr,(natom*nsppol*nspinor))
 641    iall=0
 642    do iatom=1,natom
 643      if(self%oper(1)%matlu(iatom)%lpawu.ne.-1) then
 644        ndim=2*self%oper(1)%matlu(iatom)%lpawu+1
 645        ABI_ALLOCATE(s_r,(ndim,ndim,nspinor,nspinor))
 646        ABI_ALLOCATE(s_i,(ndim,ndim,nspinor,nspinor))
 647 !       write(std_out,*) "print_self",ndim
 648        call int2char4(iatom,tag_at)
 649        ABI_CHECK((tag_at(1:1)/='#'),'Bug: string length too short!')
 650        do isppol=1,nsppol
 651          write(tag_is,'(i1)')isppol
 652 !         do ispinor=1,nspinor
 653            iall=iall+1
 654 !           write(tag_is2,'(i1)')ispinor
 655 
 656 !          ===========================
 657 !           == Create name for file
 658 !          ===========================
 659            if(self%w_type=="real") then
 660              tmpfil = trim(paw_dmft%filapp)//'Self_ra-omega_iatom'//trim(tag_at)//'_isppol'//tag_is
 661            else
 662              tmpfil = trim(paw_dmft%filapp)//'Self-omega_iatom'//trim(tag_at)//'_isppol'//tag_is
 663              if(present(opt_char)) then
 664                tmpfil = trim(paw_dmft%filapp)//'Self_ra-omega_iatom'//trim(tag_at)//'_isppol'//tag_is//opt_char
 665              endif
 666            endif
 667            if(optrw==1) write(message,'(3a)') ch10,"  == Read  self function and Fermi Level on file ",trim(tmpfil)
 668            if(optrw==2) write(message,'(3a)') ch10,"  == Write self function and Fermi Level on file ",trim(tmpfil)
 669            call wrtout(std_out,message,'COLL')
 670            !unitselffunc_arr(iall)=300+iall-1
 671            unitselffunc_arr(iall) = get_unit()
 672            ABI_CHECK(unitselffunc_arr(iall) > 0, "Cannot find free IO unit!")
 673 
 674 !           write(std_out,*) "1"
 675 
 676 !          ===========================
 677 !           == Read: check that the file exists
 678 !          ===========================
 679            if(optrw==1) then
 680 !           write(std_out,*) "3"
 681              inquire(file=trim(tmpfil),exist=lexist,recl=nrecl)
 682              if((.not.lexist)) then
 683 !           write(std_out,*) "4"
 684                iexist2=0
 685                write(message,'(4x,a,i5,3a)') "File number",unitselffunc_arr(iall),&
 686 &               " called ",trim(tmpfil)," does not exist"
 687 !               write(std_out,*) lexist,nrecl
 688                call wrtout(std_out,message,'COLL')
 689              endif
 690            endif
 691            !write(std_out,*) "2"
 692 
 693 !          ===========================
 694 !           == Open file
 695 !          ===========================
 696            if(optrw==2.or.(optrw==1.and.iexist2==1)) then
 697              !write(std_out,*) "5"
 698 #ifdef FC_NAG
 699              open (unit=unitselffunc_arr(iall),file=trim(tmpfil),status='unknown',form='formatted',recl=ABI_RECL)
 700 #else
 701              open (unit=unitselffunc_arr(iall),file=trim(tmpfil),status='unknown',form='formatted')
 702 #endif
 703              rewind(unitselffunc_arr(iall))
 704              !write(std_out,*) "61",nrecl
 705              if(prtopt>=3) then
 706                write(message,'(a,a,a,i4)') 'opened file : ', trim(tmpfil), ' unit', unitselffunc_arr(iall)
 707                call wrtout(std_out,message,'COLL')
 708              endif
 709            endif
 710            !write(std_out,*) "6",nrecl
 711 
 712 !          ===========================
 713 !           == Check Header
 714 !          ===========================
 715            if(optrw==2) then
 716              write(message,'(3a,5i5,2x,e25.17)') "# natom,nsppol,nspinor,ndim,nw,fermilevel",ch10&
 717 &             ,"####",natom,nsppol,nspinor,ndim,self%nw,paw_dmft%fermie
 718              call wrtout(unitselffunc_arr(iall),message,'COLL')
 719            else if(optrw==1.and.iexist2==1) then
 720              read(unitselffunc_arr(iall),*) 
 721              read(unitselffunc_arr(iall),*,iostat=ioerr)&
 722 &              chtemp,natom_read,nsppol_read,nspinor_read,ndim_read,nw_read,fermie_read
 723              if(ioerr<0) then
 724 !              write(std_out,*)" HEADER IOERR"
 725 !              write(std_out,'(a4,2x,31(e15.8,2x))') chtemp,natom_read,nsppol_read,nspinor_read,ndim_read,nw_read,fermie_read
 726              endif
 727              if(ioerr==0) then
 728                write(message,'(a,3x,3a,i12,2a,i11,2a,i10,2a,i13,2a,i11,2a,e25.8)') ch10,"Data in Self file corresponds to",&
 729 &               ch10,"     natom",natom_read,&
 730 &               ch10,"     nsppol",nsppol_read,&
 731 &               ch10,"     nspinor",nspinor_read,&
 732 &               ch10,"     ndim",ndim_read, &
 733 &               ch10,"     nw",nw_read, &
 734 &               ch10,"     Fermi level",fermie_read 
 735                call wrtout(std_out,message,'COLL')
 736                if((natom/=natom_read).or.(nsppol_read/=nsppol).or.&
 737 &                (nspinor/=nspinor_read).or.(nw_read/=self%nw)) then
 738                  message = "Dimensions in self are not correct"
 739                  MSG_WARNING(message)
 740                  iexist2=2
 741                endif
 742              else
 743                MSG_WARNING("Self file is empty")
 744              endif
 745            endif
 746            !write(std_out,*) "7"
 747 
 748 !          ===========================
 749 !           == Write/Read self in the file
 750 !          ===========================
 751            do ifreq=1,self%nw
 752              if(optrw==2) then
 753 !               write(std_out,'(a,2x,31(e15.8,2x))') &
 754 !&              "SETEST",self%omega(ifreq),&
 755 !&              (self%oper(ifreq)%matlu(iatom)%mat(im,im,isppol,ispinor,ispinor)&
 756 !&               ,im=1,ndim)
 757 !               write(std_out,*) self%omega(ifreq),&
 758 !&              ((self%oper(ifreq)%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor)&
 759 !&               ,im=1,ndim),im1=1,ndim)
 760 !               write(message,'(2x,393(e25.17,2x))')  self%omega(ifreq),&
 761 !&              ((((self%oper(ifreq)%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)&
 762 !&              ,im=1,ndim),im1=1,ndim),ispinor=1,nspinor),ispinor1=1,nspinor)
 763 
 764                !MGNAG Runtime Error: wrtout_cpp.f90, line 896: Buffer overflow on output
 765                !Is it possible to rewrite the code below to avoid such a long message
 766                !What about Netcdf binary files ?
 767                if(nspinor==1) then
 768                  write(message,'(2x,393(e25.17,2x))')  self%omega(ifreq),&
 769 &                ((((real(self%oper(ifreq)%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)),&
 770 &                aimag(self%oper(ifreq)%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)),&
 771 &                im=1,ndim),im1=1,ndim),ispinor=1,nspinor),ispinor1=1,nspinor)
 772                else
 773                  write(message,'(2x,393(e18.10,2x))')  self%omega(ifreq),&
 774 &                ((((real(self%oper(ifreq)%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)),&
 775 &                aimag(self%oper(ifreq)%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)),&
 776 &                im=1,ndim),im1=1,ndim),ispinor=1,nspinor),ispinor1=1,nspinor)
 777                endif
 778                call wrtout(unitselffunc_arr(iall),message,'COLL')
 779 !               write(std_out,*) unitselffunc_arr(iall)
 780              else if(optrw==1.and.iexist2==1.and.ioerr==0) then
 781            !write(std_out,*) "8"
 782 !               read(unitselffunc_arr(iall),'(2x,31(e15.8,2x))',iostat=ioerr) &
 783 !&              xtemp,(s_r(im),s_i(im),im=1,ndim)
 784                read(unitselffunc_arr(iall),*,iostat=ioerr) xtemp,&
 785 &   ((((s_r(im,im1,ispinor,ispinor1),s_i(im,im1,ispinor,ispinor1),im=1,ndim),im1=1,ndim),ispinor=1,nspinor),ispinor1=1,nspinor)
 786 !               if(ioerr<0) then
 787 !                write(std_out,*)" SELF IOERR<"
 788 !               else if(ioerr>0) then
 789 !                write(std_out,*)" SELF IOERR>"
 790 !                write(std_out,'(a4,2x,31(e15.8,2x))') xtemp,(s_r(im),s_i(im),im=1,ndim)
 791 !               endif
 792                do im=1,ndim
 793                  do im1=1,ndim
 794                    do ispinor=1,nspinor
 795                      do ispinor1=1,nspinor
 796                         self%oper(ifreq)%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)&
 797 &                       =cmplx(s_r(im,im1,ispinor,ispinor1),s_i(im,im1,ispinor,ispinor1),kind=dp)
 798                      enddo
 799                    enddo
 800                  enddo
 801                enddo
 802              endif
 803            enddo ! ifreq
 804 !          ===========================
 805 !           == Write/Read hdc in the file
 806 !          ===========================
 807            if(optrw==2) then
 808 !             write(std_out,'(a,2x,31(e15.8,2x))') &
 809 !&            "SETEST #dc ",(self%hdc%matlu(iatom)%mat(im,im,isppol,ispinor,ispinor),im=1,ndim)
 810              write(message,'(a,2x,31(e25.17,2x))') &
 811 &            "#dc ",((self%hdc%matlu(iatom)%mat(im,im,isppol,ispinor,ispinor),im=1,ndim),ispinor=1,nspinor)
 812              call wrtout(unitselffunc_arr(iall),message,'COLL')
 813            else if(optrw==1.and.iexist2==1.and.ioerr==0) then
 814          !write(std_out,*) "8"
 815              read(unitselffunc_arr(iall),*,iostat=ioerr) &
 816 &             chtemp,((s_r(im,1,ispinor,1),s_i(im,1,ispinor,1),im=1,ndim),ispinor=1,nspinor)
 817              if(ioerr<0) then
 818 !              write(std_out,*)" HDC IOERR<",ioerr
 819              else if(ioerr>0) then
 820 !              write(std_out,*)" HDC IOERR>",ioerr
 821              endif
 822              do ispinor=1,nspinor
 823                do im=1,ndim
 824                  self%hdc%matlu(iatom)%mat(im,im,isppol,ispinor,ispinor)&
 825 &                 =cmplx(s_r(im,1,ispinor,1),s_i(im,1,ispinor,1),kind=dp)
 826                enddo
 827              enddo
 828            endif
 829            close(unitselffunc_arr(iall))
 830 !         enddo ! ispinor
 831        enddo ! isppol
 832        ABI_DEALLOCATE(s_r)
 833        ABI_DEALLOCATE(s_i)
 834      endif ! lpawu=/-1
 835    enddo ! iatom
 836    ABI_DEALLOCATE(unitselffunc_arr)
 837  endif ! optrw==2.or.myproc==master
 838 ! call xmpi_barrier(spacecomm)
 839            !write(std_out,*) "9"
 840 
 841 !  ===========================
 842 !  == Error messages 
 843 !  ===========================
 844  if(optrw==1) then
 845 !   call xmpi_barrier(spacecomm)
 846    ncount=natom*nsppol*(nspinor**2)*(self%nw+1)*(maxval(self%oper(1)%matlu(:)%lpawu)*2+1)**2
 847    !write(std_out,*) ncount,maxval(pawtab(:)%lpawu)*2+1
 848    call xmpi_bcast(iexist2,master,spacecomm ,ier)
 849    call xmpi_bcast(ioerr,master,spacecomm ,ier)
 850    if(iexist2==0.or.ioerr<0.or.ioerr>0) then
 851      message = "Self file does not exist or is incomplete"
 852      MSG_WARNING(message)
 853      if(iexist2==0) then
 854        write(message,'(4x,2a)') "File does not exist"
 855        call wrtout(std_out,message,'COLL')
 856      endif
 857      if(ioerr<0) then
 858        write(message,'(4x,2a)') "End of file reached"
 859        call wrtout(std_out,message,'COLL')
 860      endif
 861      if(ioerr>0) then
 862        write(message,'(4x,2a)') "Error during read statement"
 863        call wrtout(std_out,message,'COLL')
 864      endif
 865      if(paw_dmft%dmft_solv/=4) then
 866        write(message,'(4x,a,a,5i5,2x,e14.7)') "-> Put Self-Energy Equal to double counting term"
 867      else if(paw_dmft%dmft_solv==4) then
 868        write(message,'(4x,a,a,5i5,2x,e14.7)') "-> Put Self-Energy Equal to dc term - shift"
 869        call wrtout(std_out,message,'COLL')
 870        write(message,'(4x,a,a,5i5,2x,e14.7)') " No self energy is given, change dmft_rslf"
 871        MSG_ERROR(message)
 872      endif
 873      call wrtout(std_out,message,'COLL')
 874      do ifreq=1,self%nw
 875 !       write(std_out,*) "before",self%oper(1)%matlu(1)%mat(1,1,1,1,1)
 876 !       write(std_out,*) "before",self%hdc%matlu(1)%mat(1,1,1,1,1)
 877        call copy_matlu(self%hdc%matlu,self%oper(ifreq)%matlu,natom)
 878 !       write(std_out,*) "after",self%oper(1)%matlu(1)%mat(1,1,1,1,1)
 879 !       write(std_out,*) "before",self%hdc%matlu(1)%mat(1,1,1,1,1)
 880        if(paw_dmft%dmft_solv==4) then
 881 !         if(ifreq==1) write(std_out,*) "shift",self%qmc_shift(1)
 882          call shift_matlu(self%oper(ifreq)%matlu,natom,cmplx(self%qmc_shift,0.d0,kind=dp),1)
 883 !         if(ifreq==1) write(std_out,*) "self after dc and shift",self%oper(ifreq)%matlu(1)%mat(1,1,1,1,1)
 884 !         if(ifreq==1) write(std_out,*) "shift",self%qmc_shift(1)
 885        endif
 886      enddo
 887    else ! test read successfull
 888 !   call xmpi_barrier(spacecomm)
 889 !! lignes 924-928 semblent inutiles puisque la valeur de paw_dmft%fermie creee
 890 !! en ligne 927 est ecrasee en ligne 992. BA+jmb 
 891 !!     ABI_ALLOCATE(fermie_read2,(1))
 892 !!     fermie_read2(1)=fermie_read
 893 !!     call xmpi_sum(fermie_read2,spacecomm ,ier)
 894 !!     paw_dmft%fermie=fermie_read2(1)
 895 !!     ABI_DEALLOCATE(fermie_read2)
 896   
 897 !  ===========================
 898 !   bcast to other proc
 899 !  ===========================
 900      ABI_ALLOCATE(buffer,(ncount))
 901      ABI_ALLOCATE(fermie_read2,(1))
 902      buffer(:)=czero
 903 !! BA+jmb
 904      fermie_read2=zero
 905    !write(std_out,*) self%nw
 906      if(myproc==master) then
 907      
 908 !               == Send read data to all process
 909        icount=0
 910        fermie_read2(1)=fermie_read
 911        do iatom=1,natom
 912          if(self%oper(1)%matlu(iatom)%lpawu.ne.-1) then
 913            ndim=2*self%oper(1)%matlu(iatom)%lpawu+1
 914            do isppol=1,nsppol
 915              do ispinor=1,nspinor
 916 !               Self energy-----------
 917                do ifreq=1,self%nw
 918                  do ispinor1=1,nspinor
 919                    do im=1,ndim
 920                      do im1=1,ndim
 921                        icount=icount+1
 922                        if(icount.gt.ncount) then
 923                          write(message,'(2a,2i5)') ch10,"Error buffer",icount,ncount
 924                          iexit=1
 925                          MSG_ERROR(message)
 926                        endif
 927                        buffer(icount)=self%oper(ifreq)%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)
 928                      enddo
 929                    enddo
 930                  enddo  ! ispinor1
 931                enddo
 932 !               Double counting-------
 933                do im=1,ndim
 934                  icount=icount+1
 935                  if(icount.gt.ncount) then
 936                    write(message,'(2a,2i5)') ch10,"Error buffer",icount,ncount
 937                    iexit=1
 938                    MSG_ERROR(message)
 939                  endif
 940                  buffer(icount)=self%hdc%matlu(iatom)%mat(im,im,isppol,ispinor,ispinor)
 941                enddo
 942              enddo  ! ispinor
 943            enddo ! isppol
 944          endif ! lpawu=/-1
 945        enddo ! iatom
 946      endif
 947 !    call xmpi_bcast(buffer,master,spacecomm ,ier)
 948 !    call xmpi_sum(iexit,spacecomm ,ier)
 949 !!JB call xmpi_barrier(spacecomm) 
 950      call xmpi_sum(buffer,spacecomm ,ier)
 951 !!JB call xmpi_barrier(spacecomm)
 952 
 953 ! bcast fermi level
 954    call xmpi_sum(fermie_read2,spacecomm ,ier)
 955 
 956      if(ier/=0) then
 957        message =  "error in xmpi_sum in rw_self"
 958        MSG_ERROR(message)
 959      endif
 960      paw_dmft%fermie=fermie_read2(1)
 961 !     write(std_out,*) "Fermi level",paw_dmft%fermie
 962      icount=0
 963      do iatom=1,natom
 964        if(self%oper(1)%matlu(iatom)%lpawu.ne.-1) then
 965          ndim=2*self%oper(1)%matlu(iatom)%lpawu+1
 966          do isppol=1,nsppol
 967            do ispinor=1,nspinor
 968 !             self ---------------
 969              do ifreq=1,self%nw
 970                do ispinor1=1,nspinor
 971                  do im=1,ndim
 972                    do im1=1,ndim
 973                      icount=icount+1
 974                      self%oper(ifreq)%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)=buffer(icount)
 975                    enddo
 976                  enddo
 977                enddo
 978              enddo
 979 !             hdc  ---------------
 980              do im=1,ndim
 981                icount=icount+1
 982                self%hdc%matlu(iatom)%mat(im,im,isppol,ispinor,ispinor)=buffer(icount)
 983              enddo
 984            enddo
 985          enddo ! isppol
 986        endif ! lpawu=/-1
 987      enddo ! iatom
 988      ABI_DEALLOCATE(fermie_read2)
 989      ABI_DEALLOCATE(buffer)
 990    endif  ! test read successful
 991  endif  ! optrw==1
 992 !   call flush_unit(std_out)
 993 !   MSG_ERROR("Aboring now")
 994  if(optrw==0) then
 995    if(paw_dmft%dmft_rslf==0) then
 996      if(paw_dmft%dmft_solv/=4) then
 997        write(message,'(4x,a,a,5i5,2x,e14.7)') "-> Put Self-Energy Equal to double counting term"
 998      else if(paw_dmft%dmft_solv==4) then
 999        write(message,'(4x,a,a,5i5,2x,e14.7)') "-> Put Self-Energy Equal to dc term - shift"
1000      endif
1001    else if (paw_dmft%dmft_rslf==-1) then
1002      write(message,'(4x,a,a,5i5,2x,e14.7)') "-> Put Self-Energy Equal to zero"
1003    endif
1004    call wrtout(std_out,message,'COLL')
1005    do ifreq=1,self%nw
1006      if(paw_dmft%dmft_rslf==0) then
1007        call copy_matlu(self%hdc%matlu,self%oper(ifreq)%matlu,natom)
1008       ! if(ifreq==1) write(std_out,*) "self after dc",self%oper(ifreq)%matlu(1)%mat(1,1,1,1,1)
1009        if(paw_dmft%dmft_solv==4) then
1010       !   if(ifreq==1) write(std_out,*) "shift",self%qmc_shift(1)
1011          call shift_matlu(self%oper(ifreq)%matlu,natom,cmplx(self%qmc_shift,0.d0,kind=dp),1)
1012       !   if(ifreq==1) write(std_out,*) "self after dc and shift",self%oper(ifreq)%matlu(1)%mat(1,1,1,1,1)
1013       !   if(ifreq==1) write(std_out,*) "shift",self%qmc_shift(1)
1014        endif
1015      else if (paw_dmft%dmft_rslf==-1) then
1016        do iatom=1,natom
1017          if(self%oper(1)%matlu(iatom)%lpawu.ne.-1) then
1018            ndim=2*self%oper(1)%matlu(iatom)%lpawu+1
1019            do isppol=1,nsppol
1020              do ispinor=1,nspinor
1021                do im=1,ndim
1022                  self%oper(ifreq)%matlu(iatom)%mat(im,im,isppol,ispinor,ispinor)= czero
1023                enddo
1024              enddo
1025            enddo
1026          endif
1027        enddo
1028      endif
1029    enddo
1030  endif
1031 
1032 ! write(std_out,*) "optrw,use_fixed_self,istep,iter,istep_imp,iter_imp"
1033 ! write(std_out,*) optrw,paw_dmft%use_fixed_self,istep,iter,istep_imp,iter_imp
1034  if((optrw==1.or.optrw==3).and.paw_dmft%use_fixed_self>0.and.istep<=istep_imp.and.iter<=iter_imp) then
1035    write(message,'(4x,a)') "-> Put Self-Energy Equal to imposed self-energy"
1036    call wrtout(std_out,message,'COLL')
1037    iatu=0
1038    do iatom=1,natom
1039      if(self%oper(1)%matlu(iatom)%lpawu.ne.-1) then
1040        iatu=iatu+1
1041        do ifreq=1,self%nw
1042          ndim=2*self%oper(1)%matlu(iatom)%lpawu+1
1043          do isppol=1,nsppol
1044            do ispinor=1,nspinor
1045              do im=1,ndim
1046                if(nspinor==1) then
1047                  self%oper(ifreq)%matlu(iatom)%mat(im,im,isppol,ispinor,ispinor)= paw_dmft%fixed_self(im,im,isppol,iatu)
1048 !            write(std_out,*) paw_dmft%fixed_self(im,im,isppol,iatu)
1049                else
1050                  self%oper(ifreq)%matlu(iatom)%mat(im,im,isppol,ispinor,ispinor)= paw_dmft%fixed_self(im,im,ispinor,iatu)
1051 !                 write(message,'(a,i4,i4,2x,e20.10)') " Fixed self not implemented for nspinor==2"
1052 !                 call wrtout(std_out,  message,'COLL')
1053 !                 MSG_ERROR("Aboring now")
1054                endif
1055              enddo
1056            enddo
1057          enddo
1058        enddo
1059      endif
1060    enddo
1061 
1062  endif
1063 
1064 
1065 
1066 
1067 end subroutine rw_self

m_self/self_type [ Types ]

[ Top ] [ m_self ] [ Types ]

NAME

  self_type

FUNCTION

  This structured datatype contains the necessary data

SOURCE

62  type, public :: self_type ! for each atom
63 
64   integer :: dmft_nwlo
65   ! dmft frequencies
66 
67   character(len=4) :: w_type
68   ! type of frequencies used
69 
70   integer :: nw
71   ! dmft frequencies
72 
73   integer :: iself_cv
74   ! integer for self convergence
75 
76   integer :: dmft_nwli
77   ! dmft frequencies
78 
79   real(dp), pointer :: omega(:) => null()
80   ! value of frequencies
81 
82   real(dp), allocatable :: qmc_shift(:)
83   ! value of frequencies
84 
85   real(dp), allocatable :: qmc_xmu(:)
86   ! value of frequencies
87 
88   type(oper_type), allocatable :: oper(:)
89   ! self function  in different basis
90 
91   type(oper_type):: hdc
92   ! self function  in different basis
93 
94  end type self_type