TABLE OF CONTENTS
ABINIT/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
paw_ij_free,paw_ij_init,paw_ij_nullify,pawpupot,wrtout
SOURCE
36 #if defined HAVE_CONFIG_H 37 #include "config.h" 38 #endif 39 40 41 #include "abi_common.h" 42 43 subroutine ldau_self(cryst_struc,green,paw_dmft,pawtab,self,opt_ldau,prtopt) 44 45 use m_profiling_abi 46 47 use defs_basis 48 use defs_datatypes 49 use m_errors 50 use m_crystal, only : crystal_t 51 use m_green, only : green_type 52 use m_self, only : self_type 53 use m_oper, only : print_oper 54 use m_energy, only : energy_type,compute_energy 55 56 use m_pawtab, only : pawtab_type 57 use m_paw_ij, only : paw_ij_type, paw_ij_init, paw_ij_free, paw_ij_nullify 58 use m_pawdij, only : pawpupot 59 use m_paw_dmft, only : paw_dmft_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 'ldau_self' 65 use interfaces_14_hidewrite 66 !End of the abilint section 67 68 implicit none 69 70 !Arguments ------------------------------------ 71 !scalars 72 ! type(pawang_type), intent(in) :: pawang 73 type(crystal_t),intent(in) :: cryst_struc 74 type(green_type), intent(in) :: green 75 type(paw_dmft_type), intent(in) :: paw_dmft 76 type(pawtab_type),intent(in) :: pawtab(cryst_struc%ntypat) 77 type(self_type), intent(inout) :: self !vz_i 78 integer, intent(in) :: opt_ldau,prtopt 79 80 !Local variables ------------------------------ 81 !scalars 82 character(len=500) :: message 83 integer :: iatom,idijeff,isppol,ispinor,ispinor1,itypat,lpawu 84 integer :: ifreq,im,im1,ldim,natom,nspden,nsppol,nspinor,nsploop 85 !arrays 86 integer,parameter :: spinor_idxs(2,4)=RESHAPE((/1,1,2,2,1,2,2,1/),(/2,4/)) 87 real(dp),allocatable :: vpawu(:,:,:,:) 88 type(paw_ij_type), allocatable :: paw_ij(:) 89 90 !************************************************************************ 91 92 natom=cryst_struc%natom 93 nsppol=paw_dmft%nsppol 94 nspinor=paw_dmft%nspinor 95 nspden=paw_dmft%nspden 96 nsploop=max(paw_dmft%nsppol,paw_dmft%nspinor**2) 97 isppol=0 98 ispinor=0 99 ispinor1=0 100 !write(std_out,*) "nsploop ",nsploop 101 if(opt_ldau==1) then 102 end if 103 ABI_DATATYPE_ALLOCATE(paw_ij,(cryst_struc%natom)) 104 call paw_ij_nullify(paw_ij) 105 call paw_ij_init(paw_ij,2,paw_dmft%nspinor,paw_dmft%nsppol,paw_dmft%nspden, & 106 & 1,paw_dmft%natom,cryst_struc%ntypat,cryst_struc%typat,pawtab,has_pawu_occ=1) 107 108 !=============================== 109 !Impose density matrix 110 !=============================== 111 !todo_ab build pawrhoij or suppress it in setnoccmmp 112 !dimdmat=2*maxval(pawtab(:)%lpawu)+1 113 !call setnoccmmp(0,dimdmat,paw_dmft%dmatpawu(1:dimdmat,1:dimdmat,1:nsppol*nspinor,1:natpawu),& 114 !& 1,1,cryst_struc%indsym,cryst_struc%natom,cryst_struc%natom,natpawu,& 115 !& paw_dmft%nspinor,paw_dmft%nsppol,cryst_struc%nsym,cryst_struc%ntypat,paw_ij,pawang,3,& 116 !& pawrhoij_dum,pawtab,cryst_struc%spinat,cryst_struc%symafm,struct%typat,0,10) 117 118 do iatom=1,cryst_struc%natom 119 itypat=cryst_struc%typat(iatom) 120 lpawu=paw_dmft%lpawu(iatom) 121 if(lpawu.ne.-1) then 122 ldim=2*lpawu+1 123 ABI_ALLOCATE(vpawu,(paw_ij(iatom)%cplex_dij,ldim,ldim,paw_ij(iatom)%ndij)) 124 125 paw_ij(iatom)%nocctot(:)=zero ! contains nmmp in the n m representation 126 ! =============================== 127 ! Begin loop over spin/spinors to initialize noccmmp 128 do idijeff=1,nsploop 129 ! =============================== 130 if(nsploop==2) then 131 isppol=spinor_idxs(1,idijeff) 132 ispinor=1 133 ispinor1=1 134 else if(nsploop==4) then 135 isppol=1 136 ispinor=spinor_idxs(1,idijeff) 137 ispinor1=spinor_idxs(2,idijeff) 138 else if(nsploop==1) then 139 isppol=1 140 ispinor=1 141 ispinor1=1 142 else 143 write(message,'(2a)') " BUG in ldau_self: nsploop should be equal to 1, 2 or 4" 144 call wrtout(std_out,message,'COLL') 145 end if 146 ! =============================== 147 ! Initialize noccmmp 148 ! =============================== 149 do im1 = 1 , ldim 150 do im = 1 , ldim 151 ! paw_ij(iatom)%noccmmp(1,im,im1,idijeff)=real(green%occup%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)) 152 ! paw_ij(iatom)%noccmmp(2,im,im1,idijeff)=imag(green%occup%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)) 153 paw_ij(iatom)%noccmmp(1,im,im1,idijeff)=real(green%occup%matlu(iatom)%mat(im1,im,isppol,ispinor,ispinor1)) 154 paw_ij(iatom)%noccmmp(2,im,im1,idijeff)=aimag(green%occup%matlu(iatom)%mat(im1,im,isppol,ispinor,ispinor1)) 155 end do 156 end do 157 158 ! paw_ij(1)%noccmmp(1,3,3,1)=0.d0 159 ! paw_ij(1)%noccmmp(1,5,5,1)=0.d0 160 ! paw_ij(2)%noccmmp(1,3,3,2)=0.d0 161 ! paw_ij(2)%noccmmp(1,5,5,2)=0.d0 162 ! =============================== 163 ! Compute nocctot for pawpupot = 164 ! =============================== 165 do im1=1,ldim 166 if(nsploop==4) then 167 if(idijeff<=2) then 168 paw_ij(iatom)%nocctot(1)=paw_ij(iatom)%nocctot(1)+paw_ij(iatom)%noccmmp(1,im1,im1,idijeff) 169 end if 170 else 171 paw_ij(iatom)%nocctot(idijeff)=paw_ij(iatom)%nocctot(idijeff)+paw_ij(iatom)%noccmmp(1,im1,im1,idijeff) 172 end if 173 end do 174 ! write(message,'(2a)') ch10," == The noccmmp matrix is" 175 ! call wrtout(std_out,message,'COLL') 176 ! do im=1,ldim 177 ! write(message,'(12(1x,9(1x,f10.5)))') (paw_ij(iatom)%noccmmp(1,im,im1,idijeff),im1=1,ldim) 178 ! call wrtout(std_out,message,'COLL') 179 ! write(message,'(12(1x,9(1x,f10.5)))') (paw_ij(iatom)%noccmmp(2,im,im1,idijeff),im1=1,ldim) 180 ! call wrtout(std_out,message,'COLL') 181 ! end do 182 ! ! write(message,'(2a)') ch10," == The nocctot are" 183 ! call wrtout(std_out,message,'COLL') 184 ! write(std_out,*) paw_ij(iatom)%nocctot(idijeff) 185 end do 186 paw_ij(iatom)%has_pawu_occ=2 187 188 ! warning dmft works here if nspden=nsppol (this is checked elsewhere) 189 190 ! =============================== 191 ! Compute LDA+U vpawu from noccmmp 192 ! =============================== 193 call pawpupot(paw_ij(iatom)%cplex_dij,paw_ij(iatom)%ndij,& 194 & paw_ij(iatom)%noccmmp,paw_ij(iatom)%nocctot,& 195 & 2,pawtab(itypat),vpawu) 196 ! do idijeff=1,size(vpawu,4) 197 ! write(message,'(2a)') ch10," == The vpawu matrix is" 198 ! call wrtout(std_out,message,'COLL') 199 ! do im=1,ldim 200 ! write(message,'(12(1x,9(1x,f10.5)))') (vpawu(1,im,im1,idijeff),im1=1,ldim) 201 ! call wrtout(std_out,message,'COLL') 202 ! write(message,'(12(1x,9(1x,f10.5)))') (vpawu(2,im,im1,idijeff),im1=1,ldim) 203 ! call wrtout(std_out,message,'COLL') 204 ! end do 205 ! end do 206 207 ! =============================== 208 ! Begin loop over spin/spinors to compute self%oper 209 do idijeff=1,nsploop 210 ! =============================== 211 if(nsploop==2) then 212 isppol=spinor_idxs(1,idijeff) 213 ispinor=1 214 ispinor1=1 215 else if(nsploop==4) then 216 isppol=1 217 ispinor=spinor_idxs(1,idijeff) 218 ispinor1=spinor_idxs(2,idijeff) 219 else if(nsploop==1) then 220 isppol=1 221 ispinor=1 222 ispinor1=1 223 else 224 write(message,'(2a)') " BUG in ldau_self: nsploop should be equal to 1, 2 or 4" 225 call wrtout(std_out,message,'COLL') 226 end if 227 228 ! =============================== 229 ! vpawu -> self%oper 230 ! =============================== 231 do im1 = 1 , ldim 232 do im = 1 , ldim 233 do ifreq=1,self%nw 234 self%oper(ifreq)%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)=& 235 ! & cmplx(vpawu(1,im1,im),vpawu(2,im1,im),kind=dp) 236 ! One take the transpose in orbital index to be coherent with the 237 ! current DFT+U implementation in Abinit. 238 & cmplx(vpawu(1,im,im1,idijeff),vpawu(2,im,im1,idijeff),kind=dp) 239 end do 240 end do 241 end do 242 243 end do ! idijeff 244 ! write(std_out,*) "check im,im1 in vpawu",iatom 245 ABI_DEALLOCATE(vpawu) 246 ! =============================== 247 ! Compute energy 248 ! =============================== 249 250 end if 251 end do 252 253 call paw_ij_free(paw_ij) 254 ABI_DATATYPE_DEALLOCATE(paw_ij) 255 256 if(abs(prtopt)>0) then 257 end if 258 259 end subroutine ldau_self