Branch data 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 : : !!***
|
