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