TABLE OF CONTENTS


ABINIT/m_ldau_self [ Modules ]

[ Top ] [ Modules ]

NAME

  m_ldau_self

FUNCTION

 Compute DFT+U self energy for DMFT

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 #include "abi_common.h"
29 
30 MODULE m_ldau_self
31 
32  use defs_basis
33 
34  implicit none
35 
36  private
37 
38  public :: ldau_self

m_ldau_self/ldau_self [ Functions ]

[ Top ] [ m_ldau_self ] [ Functions ]

NAME

 ldau_self

FUNCTION

 Use LDA+U to compute self-energy

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

  cryst_struc
  istep    =  step of iteration for LDA.
  mpi_enreg=informations about MPI parallelization
  paw_dmft =  data for self-consistent LDA+DMFT calculations.
  pawang <type(pawang)>=paw angular mesh and related data
  pawtab <type(pawtab)>

OUTPUT

NOTES

PARENTS

      impurity_solve,spectral_function

CHILDREN

      pawpupot,wrtout

SOURCE

 77 subroutine ldau_self(cryst_struc,green,paw_dmft,pawtab,self,opt_ldau,prtopt)
 78 
 79  use m_abicore
 80 
 81  use defs_basis
 82  use defs_datatypes
 83  use m_errors
 84  use m_crystal, only : crystal_t
 85  use m_green, only : green_type
 86  use m_self, only : self_type
 87  use m_oper, only : print_oper
 88  use m_energy, only : energy_type,compute_energy
 89 
 90  use m_pawtab, only : pawtab_type
 91  use m_pawdij, only : pawpupot
 92  use m_paw_dmft, only : paw_dmft_type
 93  use m_paw_correlations, only : setnoccmmp
 94 
 95 !This section has been created automatically by the script Abilint (TD).
 96 !Do not modify the following lines by hand.
 97 #undef ABI_FUNC
 98 #define ABI_FUNC 'ldau_self'
 99 !End of the abilint section
100 
101  implicit none
102 
103 !Arguments ------------------------------------
104 !scalars
105 ! type(pawang_type), intent(in) :: pawang
106  type(crystal_t),intent(in) :: cryst_struc
107  type(green_type), intent(in) :: green
108  type(paw_dmft_type), intent(in)  :: paw_dmft
109  type(pawtab_type),target,intent(in)  :: pawtab(cryst_struc%ntypat)
110  type(self_type), intent(inout) :: self !vz_i
111  integer, intent(in) :: opt_ldau,prtopt
112 
113 !Local variables ------------------------------
114 !scalars
115  character(len=500) :: message
116  integer :: iatom,idijeff,isppol,ispinor,ispinor1,lpawu
117  integer :: ifreq,im,im1,ldim,natom,nocc,nsppol,nspinor,nsploop
118 !arrays
119  integer,parameter :: spinor_idxs(2,4)=RESHAPE((/1,1,2,2,1,2,2,1/),(/2,4/))
120  real(dp),allocatable :: noccmmp(:,:,:,:),nocctot(:)
121  real(dp),allocatable :: vpawu(:,:,:,:)
122  type(pawtab_type),pointer :: pawtab_
123 
124 !************************************************************************
125 
126  natom=cryst_struc%natom
127  nsppol=paw_dmft%nsppol
128  nspinor=paw_dmft%nspinor
129  nsploop=max(nsppol,nspinor**2)
130  nocc=nsploop
131  isppol=0
132  ispinor=0
133  ispinor1=0
134 
135 !===============================
136 !Impose density matrix
137 !===============================
138 !todo_ab build pawrhoij or suppress it in setnoccmmp
139 !dimdmat=2*maxval(pawtab(:)%lpawu)+1
140 !call setnoccmmp(0,dimdmat,paw_dmft%dmatpawu(1:dimdmat,1:dimdmat,1:nsppol*nspinor,1:natpawu),&
141 !&   1,1,cryst_struc%indsym,cryst_struc%natom,cryst_struc%natom,natpawu,&
142 !&   paw_dmft%nspinor,paw_dmft%nsppol,cryst_struc%nsym,cryst_struc%ntypat,paw_ij,pawang,3,&
143 !&   pawrhoij_dum,pawtab,cryst_struc%spinat,cryst_struc%symafm,struct%typat,0,10)
144 
145  do iatom=1,cryst_struc%natom
146    lpawu=paw_dmft%lpawu(iatom)
147    pawtab_ => pawtab(cryst_struc%typat(iatom))
148    if(lpawu.ne.-1) then
149      ldim=2*lpawu+1
150      ABI_ALLOCATE(vpawu,(2,ldim,ldim,nocc))
151 
152      ABI_ALLOCATE(noccmmp,(2,2*pawtab_%lpawu+1,2*pawtab_%lpawu+1,nocc))
153      ABI_ALLOCATE(nocctot,(nocc))
154      noccmmp(:,:,:,:)=zero ; nocctot(:)=zero ! contains nmmp in the n m representation
155 
156 !    ===============================
157 !    Begin loop over spin/spinors to initialize noccmmp
158      do idijeff=1,nsploop
159 !      ===============================
160        if(nsploop==2) then
161          isppol=spinor_idxs(1,idijeff)
162          ispinor=1
163          ispinor1=1
164        else if(nsploop==4) then
165          isppol=1
166          ispinor=spinor_idxs(1,idijeff)
167          ispinor1=spinor_idxs(2,idijeff)
168        else if(nsploop==1) then
169          isppol=1
170          ispinor=1
171          ispinor1=1
172        else
173          write(message,'(2a)') " BUG in ldau_self: nsploop should be equal to 1, 2 or 4"
174          call wrtout(std_out,message,'COLL')
175        end if
176 !      ===============================
177 !      Initialize noccmmp
178 !      ===============================
179        do im1 = 1 , ldim
180          do im = 1 ,  ldim
181 !          noccmmp(1,im,im1,idijeff)=real(green%occup%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1))
182 !          noccmmp(2,im,im1,idijeff)=imag(green%occup%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1))
183            noccmmp(1,im,im1,idijeff)=real(green%occup%matlu(iatom)%mat(im1,im,isppol,ispinor,ispinor1))
184            noccmmp(2,im,im1,idijeff)=aimag(green%occup%matlu(iatom)%mat(im1,im,isppol,ispinor,ispinor1))
185          end do
186        end do
187 
188 !      ===============================
189 !      Compute nocctot for pawpupot =
190 !      ===============================
191        do im1=1,ldim
192          if(nsploop==4) then
193            if(idijeff<=2) then
194              nocctot(1)=nocctot(1)+noccmmp(1,im1,im1,idijeff)
195            end if
196          else
197            nocctot(idijeff)=nocctot(idijeff)+noccmmp(1,im1,im1,idijeff)
198          end if
199        end do
200 !      write(message,'(2a)') ch10," == The noccmmp matrix is"
201 !      call wrtout(std_out,message,'COLL')
202 !      do im=1,ldim
203 !      write(message,'(12(1x,9(1x,f10.5)))') (noccmmp(1,im,im1,idijeff),im1=1,ldim)
204 !      call wrtout(std_out,message,'COLL')
205 !      write(message,'(12(1x,9(1x,f10.5)))') (noccmmp(2,im,im1,idijeff),im1=1,ldim)
206 !      call wrtout(std_out,message,'COLL')
207 !      end do
208 !      !     write(message,'(2a)') ch10," == The nocctot are"
209 !      call wrtout(std_out,message,'COLL')
210 !      write(std_out,*) nocctot(idijeff)
211      end do
212 
213 !    warning  dmft works here if nspden=nsppol (this is checked elsewhere)
214 
215 !    ===============================
216 !    Compute LDA+U vpawu from noccmmp
217 !    ===============================
218      call pawpupot(2,nocc,noccmmp,nocctot,2,pawtab_,vpawu)
219 !    do idijeff=1,size(vpawu,4)
220 !    write(message,'(2a)') ch10," == The vpawu matrix is"
221 !    call wrtout(std_out,message,'COLL')
222 !    do im=1,ldim
223 !    write(message,'(12(1x,9(1x,f10.5)))') (vpawu(1,im,im1,idijeff),im1=1,ldim)
224 !    call wrtout(std_out,message,'COLL')
225 !    write(message,'(12(1x,9(1x,f10.5)))') (vpawu(2,im,im1,idijeff),im1=1,ldim)
226 !    call wrtout(std_out,message,'COLL')
227 !    end do
228 !    end do
229 
230 !    ===============================
231 !    Begin loop over spin/spinors to compute self%oper
232      do idijeff=1,nsploop
233 !      ===============================
234        if(nsploop==2) then
235          isppol=spinor_idxs(1,idijeff)
236          ispinor=1
237          ispinor1=1
238        else if(nsploop==4) then
239          isppol=1
240          ispinor=spinor_idxs(1,idijeff)
241          ispinor1=spinor_idxs(2,idijeff)
242        else if(nsploop==1) then
243          isppol=1
244          ispinor=1
245          ispinor1=1
246        else
247          write(message,'(2a)') " BUG in ldau_self: nsploop should be equal to 1, 2 or 4"
248          call wrtout(std_out,message,'COLL')
249        end if
250 
251 !      ===============================
252 !      vpawu -> self%oper
253 !      ===============================
254        do im1 = 1 , ldim
255          do im = 1 ,  ldim
256            do ifreq=1,self%nw
257              self%oper(ifreq)%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)=&
258 !            &             cmplx(vpawu(1,im1,im),vpawu(2,im1,im),kind=dp)
259 !            One take the transpose in orbital index to be coherent with the
260 !            current DFT+U implementation in Abinit.
261 &             cmplx(vpawu(1,im,im1,idijeff),vpawu(2,im,im1,idijeff),kind=dp)
262            end do
263          end do
264        end do
265 
266      end do ! idijeff
267 !    write(std_out,*) "check im,im1 in vpawu",iatom
268      ABI_DEALLOCATE(vpawu)
269 
270 !    ===============================
271 !    Compute energy
272 !    ===============================
273 
274      ABI_DEALLOCATE(noccmmp)
275      ABI_DEALLOCATE(nocctot)
276    end if
277  end do
278 
279  if(abs(prtopt)>0) then
280  end if
281 
282 end subroutine ldau_self