Line data Source code
1 : !!****m* ABINIT/m_dftu_self
2 : !! NAME
3 : !! m_dftu_self
4 : !!
5 : !! FUNCTION
6 : !! Compute DFT+U self energy for DMFT
7 : !!
8 : !! COPYRIGHT
9 : !! Copyright (C) 2006-2026 ABINIT group (BAmadon)
10 : !! This file is distributed under the terms of the
11 : !! GNU General Public License, see ~abinit/COPYING
12 : !! or http://www.gnu.org/copyleft/gpl.txt .
13 : !!
14 : !! INPUTS
15 : !!
16 : !! OUTPUT
17 : !!
18 : !! SOURCE
19 :
20 : #if defined HAVE_CONFIG_H
21 : #include "config.h"
22 : #endif
23 :
24 : #include "abi_common.h"
25 :
26 : MODULE m_dftu_self
27 :
28 : use defs_basis
29 : use m_abicore
30 : use m_errors
31 :
32 : use m_crystal, only : crystal_t
33 : use m_green, only : green_type
34 : use m_self, only : self_type
35 : use m_oper, only : print_oper
36 : use m_energy, only : energy_type,compute_energy
37 : use m_pawtab, only : pawtab_type
38 : use m_pawdij, only : pawpupot
39 : use m_paw_dmft, only : paw_dmft_type
40 : use m_paw_correlations, only : setnoccmmp
41 :
42 : implicit none
43 :
44 : private
45 :
46 : public :: dftu_self
47 : !!***
48 :
49 : contains
50 :
51 : !!****f* m_dftu_self/dftu_self
52 : !! NAME
53 : !! dftu_self
54 : !!
55 : !! FUNCTION
56 : !! Use DFT+U to compute self-energy
57 : !!
58 : !! COPYRIGHT
59 : !! Copyright (C) 1999-2026 ABINIT group (BAmadon)
60 : !! This file is distributed under the terms of the
61 : !! GNU General Public License, see ~abinit/COPYING
62 : !! or http://www.gnu.org/copyleft/gpl.txt .
63 : !! For the initials of contributors, see ~abinit/doc/developers/contributors.txt .
64 : !!
65 : !! INPUTS
66 : !! cryst_struc
67 : !! istep = step of iteration for DFT.
68 : !! mpi_enreg=information about MPI parallelization
69 : !! paw_dmft = data for self-consistent DFT+DMFT calculations.
70 : !! pawang <type(pawang)>=paw angular mesh and related data
71 : !! pawtab <type(pawtab)>
72 : !!
73 : !! OUTPUT
74 : !!
75 : !! NOTES
76 : !!
77 : !! SOURCE
78 :
79 2 : subroutine dftu_self(cryst_struc,green,paw_dmft,pawtab,self,opt_dftu,prtopt)
80 :
81 : !Arguments ------------------------------------
82 : !scalars
83 : ! type(pawang_type), intent(in) :: pawang
84 : type(crystal_t),intent(in) :: cryst_struc
85 : type(green_type), intent(in) :: green
86 : type(paw_dmft_type), intent(in) :: paw_dmft
87 : type(pawtab_type),target,intent(in) :: pawtab(cryst_struc%ntypat)
88 : type(self_type), intent(inout) :: self !vz_i
89 : integer, intent(in) :: opt_dftu,prtopt
90 :
91 : !Local variables ------------------------------
92 : !scalars
93 : character(len=500) :: message
94 : integer :: iatom,idijeff,isppol,ispinor,ispinor1,lpawu
95 : integer :: ifreq,im,im1,ldim,natom,nocc,nsppol,nspinor,nsploop
96 : !arrays
97 : integer,parameter :: spinor_idxs(2,4)=RESHAPE((/1,1,2,2,1,2,2,1/),(/2,4/))
98 2 : real(dp),allocatable :: noccmmp(:,:,:,:),nocctot(:)
99 2 : real(dp),allocatable :: vpawu(:,:,:,:)
100 : type(pawtab_type),pointer :: pawtab_
101 :
102 : !************************************************************************
103 :
104 : ABI_UNUSED(opt_dftu)
105 :
106 2 : natom=cryst_struc%natom
107 2 : nsppol=paw_dmft%nsppol
108 2 : nspinor=paw_dmft%nspinor
109 2 : nsploop=max(nsppol,nspinor**2)
110 2 : nocc=nsploop
111 2 : isppol=0
112 2 : ispinor=0
113 2 : ispinor1=0
114 :
115 : !===============================
116 : !Impose density matrix
117 : !===============================
118 : !todo_ab build pawrhoij or suppress it in setnoccmmp
119 : !dimdmat=2*maxval(pawtab(:)%lpawu)+1
120 : !call setnoccmmp(0,dimdmat,paw_dmft%dmatpawu(1:dimdmat,1:dimdmat,1:nsppol*nspinor,1:natpawu),&
121 : !& 1,1,cryst_struc%indsym,cryst_struc%natom,cryst_struc%natom,natpawu,&
122 : !& paw_dmft%nspinor,paw_dmft%nsppol,cryst_struc%nsym,cryst_struc%ntypat,paw_ij,pawang,3,&
123 : !& pawrhoij_dum,pawtab,cryst_struc%spinat,cryst_struc%symafm,struct%typat,0,10)
124 :
125 6 : do iatom=1,cryst_struc%natom
126 4 : lpawu=paw_dmft%lpawu(iatom)
127 4 : pawtab_ => pawtab(cryst_struc%typat(iatom))
128 6 : if(lpawu.ne.-1) then
129 2 : ldim=2*lpawu+1
130 10 : ABI_MALLOC(vpawu,(2,ldim,ldim,nocc))
131 :
132 10 : ABI_MALLOC(noccmmp,(2,2*pawtab_%lpawu+1,2*pawtab_%lpawu+1,nocc))
133 6 : ABI_MALLOC(nocctot,(nocc))
134 332 : noccmmp(:,:,:,:)=zero ; nocctot(:)=zero ! contains nmmp in the n m representation
135 :
136 : ! ===============================
137 : ! Begin loop over spin/spinors to initialize noccmmp
138 6 : do idijeff=1,nsploop
139 : ! ===============================
140 4 : if(nsploop==2) then
141 4 : isppol=spinor_idxs(1,idijeff)
142 4 : ispinor=1
143 4 : ispinor1=1
144 0 : else if(nsploop==4) then
145 0 : isppol=1
146 0 : ispinor=spinor_idxs(1,idijeff)
147 0 : ispinor1=spinor_idxs(2,idijeff)
148 0 : else if(nsploop==1) then
149 : isppol=1
150 : ispinor=1
151 : ispinor1=1
152 : else
153 0 : write(message,'(2a)') " BUG in dftu_self: nsploop should be equal to 1, 2 or 4"
154 0 : call wrtout(std_out,message,'COLL')
155 : end if
156 : ! ===============================
157 : ! Initialize noccmmp
158 : ! ===============================
159 24 : do im1 = 1 , ldim
160 124 : do im = 1 , ldim
161 : ! noccmmp(1,im,im1,idijeff)=real(green%occup%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1))
162 : ! noccmmp(2,im,im1,idijeff)=imag(green%occup%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1))
163 100 : noccmmp(1,im,im1,idijeff)=real(green%occup%matlu(iatom)%mat(im1+(ispinor-1)*ldim,im+(ispinor1-1)*ldim,isppol))
164 120 : noccmmp(2,im,im1,idijeff)=aimag(green%occup%matlu(iatom)%mat(im1+(ispinor-1)*ldim,im+(ispinor1-1)*ldim,isppol))
165 : end do
166 : end do
167 :
168 : ! ===============================
169 : ! Compute nocctot for pawpupot =
170 : ! ===============================
171 26 : do im1=1,ldim
172 24 : if(nsploop==4) then
173 0 : if(idijeff<=2) then
174 0 : nocctot(1)=nocctot(1)+noccmmp(1,im1,im1,idijeff)
175 : end if
176 : else
177 20 : nocctot(idijeff)=nocctot(idijeff)+noccmmp(1,im1,im1,idijeff)
178 : end if
179 : end do
180 : ! write(message,'(2a)') ch10," == The noccmmp matrix is"
181 : ! call wrtout(std_out,message,'COLL')
182 : ! do im=1,ldim
183 : ! write(message,'(12(1x,9(1x,f10.5)))') (noccmmp(1,im,im1,idijeff),im1=1,ldim)
184 : ! call wrtout(std_out,message,'COLL')
185 : ! write(message,'(12(1x,9(1x,f10.5)))') (noccmmp(2,im,im1,idijeff),im1=1,ldim)
186 : ! call wrtout(std_out,message,'COLL')
187 : ! end do
188 : ! ! write(message,'(2a)') ch10," == The nocctot are"
189 : ! call wrtout(std_out,message,'COLL')
190 : ! write(std_out,*) nocctot(idijeff)
191 : end do
192 :
193 : ! warning dmft works here if nspden=nsppol (this is checked elsewhere)
194 :
195 : ! ===============================
196 : ! Compute DFT+U vpawu from noccmmp
197 : ! ===============================
198 2 : call pawpupot(2,nocc,noccmmp,nocctot,2,pawtab_,vpawu)
199 : ! do idijeff=1,size(vpawu,4)
200 : ! write(message,'(2a)') ch10," == The vpawu matrix is"
201 : ! call wrtout(std_out,message,'COLL')
202 : ! do im=1,ldim
203 : ! write(message,'(12(1x,9(1x,f10.5)))') (vpawu(1,im,im1,idijeff),im1=1,ldim)
204 : ! call wrtout(std_out,message,'COLL')
205 : ! write(message,'(12(1x,9(1x,f10.5)))') (vpawu(2,im,im1,idijeff),im1=1,ldim)
206 : ! call wrtout(std_out,message,'COLL')
207 : ! end do
208 : ! end do
209 :
210 : ! ===============================
211 : ! Begin loop over spin/spinors to compute self%oper
212 6 : do idijeff=1,nsploop
213 : ! ===============================
214 4 : if(nsploop==2) then
215 4 : isppol=spinor_idxs(1,idijeff)
216 4 : ispinor=1
217 4 : ispinor1=1
218 0 : else if(nsploop==4) then
219 0 : isppol=1
220 0 : ispinor=spinor_idxs(1,idijeff)
221 0 : ispinor1=spinor_idxs(2,idijeff)
222 0 : else if(nsploop==1) then
223 : isppol=1
224 : ispinor=1
225 : ispinor1=1
226 : else
227 0 : write(message,'(2a)') " BUG in dftu_self: nsploop should be equal to 1, 2 or 4"
228 0 : call wrtout(std_out,message,'COLL')
229 : end if
230 :
231 : ! ===============================
232 : ! vpawu -> self%oper
233 : ! ===============================
234 26 : do im1 = 1 , ldim
235 124 : do im = 1 , ldim
236 3320 : do ifreq=1,self%nw
237 : self%oper(ifreq)%matlu(iatom)%mat(im+(ispinor-1)*ldim,im1+(ispinor1-1)*ldim,isppol)=&
238 : ! & cmplx(vpawu(1,im1,im),vpawu(2,im1,im),kind=dp)
239 : ! One take the transpose in orbital index to be coherent with the
240 : ! current DFT+U implementation in Abinit.
241 3300 : & cmplx(vpawu(1,im,im1,idijeff),vpawu(2,im,im1,idijeff),kind=dp)
242 : end do
243 : end do
244 : end do
245 :
246 : end do ! idijeff
247 : ! write(std_out,*) "check im,im1 in vpawu",iatom
248 2 : ABI_FREE(vpawu)
249 :
250 : ! ===============================
251 : ! Compute energy
252 : ! ===============================
253 :
254 2 : ABI_FREE(noccmmp)
255 2 : ABI_FREE(nocctot)
256 : end if
257 : end do
258 :
259 : if(abs(prtopt)>0) then
260 : end if
261 :
262 2 : end subroutine dftu_self
263 : !!***
264 :
265 : END MODULE m_dftu_self
266 : !!***
|
