Line data Source code
1 : !--------------------------------------------------------------------------------
2 : ! Copyright (c) 2016 Peter Grünberg Institut, Forschungszentrum Jülich, Germany
3 : ! This file is part of FLEUR and available as free software under the conditions
4 : ! of the MIT license as expressed in the LICENSE file in more detail.
5 : !--------------------------------------------------------------------------------
6 :
7 : MODULE m_eparas
8 : !***********************************************************************
9 : ! Calculates qlo, enerlo and sqlo, which are needed to determine the
10 : ! new energy parameters.
11 : ! Philipp Kurz 99/04
12 : !***********************************************************************
13 : ! also the 'normal' energy parameters are now included...
14 : !
15 : ! if (l_mcd) then mcd contains mcd spectrum: first index = polarization
16 : ! second = core level ; third = band index gb.2001
17 : ! corrected to work also for multiple LO's of same l at the same atom
18 : ! gb.2005
19 : !*************** ABBREVIATIONS *****************************************
20 : ! qlo : charge density of one local orbital at the current k-point
21 : ! sqlo : qlo integrated over the Brillouin zone
22 : ! enerlo : qlo*energy integrated over the Brillouin zone
23 : !***********************************************************************
24 : !
25 : CONTAINS
26 8018 : SUBROUTINE eparas(jsp,atoms,banddos,noccbd,ev_list,fmpi,ikpt,ne,we,eig,skip_t,l_evp,eigVecCoeffs,&
27 : usdus,regCharges,dos,mcd)
28 : USE m_types
29 : use m_types_dos
30 : use m_types_mcd
31 : IMPLICIT NONE
32 : TYPE(t_usdus), INTENT(IN) :: usdus
33 : TYPE(t_mpi), INTENT(IN) :: fmpi
34 : TYPE(t_atoms), INTENT(IN) :: atoms
35 : TYPE(t_banddos), INTENT(IN) :: banddos
36 : TYPE(t_eigVecCoeffs), INTENT(IN) :: eigVecCoeffs
37 : TYPE(t_regionCharges), INTENT(INOUT) :: regCharges
38 : TYPE(t_dos), INTENT(INOUT) :: dos
39 : TYPE(t_mcd), OPTIONAL, INTENT(INOUT) :: mcd
40 : ! ..
41 : ! .. Scalar Arguments ..
42 : INTEGER, INTENT (IN) :: noccbd,jsp
43 : INTEGER, INTENT (IN) :: ne,ikpt ,skip_t
44 : LOGICAL, INTENT (IN) :: l_evp
45 : INTEGER, INTENT (IN) :: ev_list(noccbd)
46 : ! ..
47 : ! .. Array Arguments ..
48 : REAL, INTENT (IN) :: eig(:)!(input%neig),
49 : REAL, INTENT (IN) :: we(noccbd)
50 :
51 : ! ..
52 : ! .. Local Scalars ..
53 : INTEGER i,l,lo,lop ,natom,nn,ntyp,m,n_dos
54 : INTEGER nt1,nt2,lm,n,ll1,ipol,icore,index
55 : REAL fac
56 : COMPLEX suma,sumb,sumab,sumba,sumaa,sumbb
57 : ! ..
58 : ! .. Local Arrays ..
59 8018 : REAL qlo(noccbd,atoms%nlod,atoms%nlod,atoms%ntype)
60 8018 : REAL qaclo(noccbd,atoms%nlod,atoms%ntype),qbclo(noccbd,atoms%nlod,atoms%ntype)
61 8018 : LOGICAL atomTypeCovered(atoms%ntype)
62 : ! ..
63 : !
64 : !---> initialize ener, sqal, enerlo and sqlo on first call
65 : !
66 :
67 8018 : IF ((ikpt.LE.fmpi%isize).AND..NOT.l_evp) THEN
68 2004 : regCharges%ener(:,:,jsp) = 0.0
69 2004 : regCharges%sqal(:,:,jsp) = 0.0
70 1052 : regCharges%enerlo(:,:,jsp) = 0.0
71 1052 : regCharges%sqlo(:,:,jsp) = 0.0
72 13868 : dos%qal(:,:,:,ikpt,jsp) = 0.0
73 : END IF
74 :
75 25290 : atomTypeCovered(:) = .FALSE.
76 : !---> l-decomposed density for each occupied state
77 : !
78 : ! DO 140 i = (skip_t+1),ne ! this I need for all states
79 113628 : DO i = 1,ne ! skip in next loop
80 121596 : DO n_dos = 1,size(banddos%dos_typelist)
81 7968 : n=banddos%dos_typelist(n_dos)
82 7968 : atomTypeCovered(n) = .TRUE.
83 7968 : fac = 1./atoms%neq(n)
84 7968 : nt1 = atoms%firstAtom(n)
85 7968 : nt2 = nt1 + atoms%neq(n) - 1
86 39840 : DO l = 0,3
87 31872 : suma = CMPLX(0.,0.)
88 31872 : sumb = CMPLX(0.,0.)
89 31872 : ll1 = l* (l+1)
90 159360 : DO m = -l,l
91 127488 : lm = ll1 + m
92 159360 : IF (.NOT.banddos%l_mcd) THEN
93 377856 : DO natom = nt1,nt2
94 250368 : suma = suma + eigVecCoeffs%abcof(i,lm,0,natom,jsp)*CONJG(eigVecCoeffs%abcof(i,lm,0,natom,jsp))
95 377856 : sumb = sumb + eigVecCoeffs%abcof(i,lm,1,natom,jsp)*CONJG(eigVecCoeffs%abcof(i,lm,1,natom,jsp))
96 : ENDDO
97 : ELSE
98 : sumaa = CMPLX(0.,0.) ; sumab = CMPLX(0.,0.)
99 : sumbb = CMPLX(0.,0.) ; sumba = CMPLX(0.,0.)
100 0 : DO natom = nt1,nt2
101 0 : sumaa = sumaa + eigVecCoeffs%abcof(i,lm,0,natom,jsp)*CONJG(eigVecCoeffs%abcof(i,lm,0,natom,jsp))
102 0 : sumbb = sumbb + eigVecCoeffs%abcof(i,lm,1,natom,jsp)*CONJG(eigVecCoeffs%abcof(i,lm,1,natom,jsp))
103 0 : sumab = sumab + eigVecCoeffs%abcof(i,lm,0,natom,jsp) *CONJG(eigVecCoeffs%abcof(i,lm,1,natom,jsp))
104 0 : sumba = sumba + eigVecCoeffs%abcof(i,lm,1,natom,jsp) *CONJG(eigVecCoeffs%abcof(i,lm,0,natom,jsp))
105 : ENDDO
106 0 : DO icore = 1, mcd%ncore(n)
107 0 : DO ipol = 1, 3
108 0 : index = 3*(n_dos-1) + ipol
109 : mcd%mcd(index,icore,ev_list(i),ikpt,jsp) = mcd%mcd(index,icore,ev_list(i),ikpt,jsp) + fac * &
110 : (sumaa * CONJG(mcd%m_mcd(icore,lm+1,index,1))*mcd%m_mcd(icore,lm+1,index,1) + &
111 : sumbb * CONJG(mcd%m_mcd(icore,lm+1,index,2))*mcd%m_mcd(icore,lm+1,index,2) + &
112 : sumab * CONJG(mcd%m_mcd(icore,lm+1,index,2))*mcd%m_mcd(icore,lm+1,index,1) + &
113 0 : sumba * CONJG(mcd%m_mcd(icore,lm+1,index,1))*mcd%m_mcd(icore,lm+1,index,2))
114 : ENDDO
115 : ENDDO
116 : ENDIF ! end MCD
117 : ENDDO
118 31872 : dos%qal(l,n_dos,ev_list(i),ikpt,jsp) = (suma+sumb*usdus%ddn(l,n,jsp))/atoms%neq(n)
119 39840 : dos%qTot(ev_list(i),ikpt,jsp) = dos%qTot(ev_list(i),ikpt,jsp) + (suma+sumb*usdus%ddn(l,n,jsp))
120 : ENDDO
121 :
122 162154 : DO l = 4, atoms%lmax(n)
123 48576 : suma = CMPLX(0.,0.)
124 48576 : sumb = CMPLX(0.,0.)
125 48576 : ll1 = l* (l+1)
126 758496 : DO m = -l,l
127 709920 : lm = ll1 + m
128 1967616 : DO natom = nt1,nt2
129 1209120 : suma = suma + eigVecCoeffs%abcof(i,lm,0,natom,jsp)*CONJG(eigVecCoeffs%abcof(i,lm,0,natom,jsp))
130 1919040 : sumb = sumb + eigVecCoeffs%abcof(i,lm,1,natom,jsp)*CONJG(eigVecCoeffs%abcof(i,lm,1,natom,jsp))
131 : ENDDO
132 : ENDDO
133 56544 : dos%qTot(ev_list(i),ikpt,jsp) = dos%qTot(ev_list(i),ikpt,jsp) + (suma+sumb*usdus%ddn(l,n,jsp))
134 : ENDDO
135 :
136 :
137 : ENDDO
138 : ENDDO
139 :
140 25290 : DO n = 1, atoms%ntype
141 17272 : IF(atomTypeCovered(n)) THEN
142 : CYCLE
143 : END IF
144 17076 : nt1 = atoms%firstAtom(n)
145 17076 : nt2 = nt1 + atoms%neq(n) - 1
146 245056 : DO i = 1,ne ! skip in next loop
147 2244064 : DO l = 0, atoms%lmax(n)
148 2006830 : suma = CMPLX(0.,0.)
149 2006830 : sumb = CMPLX(0.,0.)
150 2006830 : ll1 = l* (l+1)
151 20756840 : DO m = -l,l
152 18750010 : lm = ll1 + m
153 39854834 : DO natom = nt1,nt2
154 19097994 : suma = suma + eigVecCoeffs%abcof(i,lm,0,natom,jsp)*CONJG(eigVecCoeffs%abcof(i,lm,0,natom,jsp))
155 37848004 : sumb = sumb + eigVecCoeffs%abcof(i,lm,1,natom,jsp)*CONJG(eigVecCoeffs%abcof(i,lm,1,natom,jsp))
156 : ENDDO
157 : ENDDO
158 2226792 : dos%qTot(ev_list(i),ikpt,jsp) = dos%qTot(ev_list(i),ikpt,jsp) + (suma+sumb*usdus%ddn(l,n,jsp))
159 : ENDDO
160 : ENDDO
161 : ENDDO
162 :
163 :
164 : !
165 : !---> perform Brillouin zone integration and summation over the
166 : !---> bands in order to determine the energy parameters for each
167 : !---> atom and angular momentum
168 : !
169 40090 : DO l = 0,3
170 109178 : DO n = 1,atoms%ntype
171 101160 : if (banddos%map_atomtype(n)==0) then
172 947640 : DO i = (skip_t+1),noccbd
173 879336 : suma = CMPLX(0.,0.)
174 879336 : sumb = CMPLX(0.,0.)
175 1776016 : DO natom = atoms%firstAtom(n), atoms%firstAtom(n) + atoms%neq(n) - 1
176 4483400 : suma=suma+dot_product(eigVecCoeffs%abcof(i,l* (l+1)-l:l* (l+1)+l,0,natom,jsp),eigVecCoeffs%abcof(i,l* (l+1)-l:l* (l+1)+l,0,natom,jsp))
177 5362736 : sumb=sumb+dot_product(eigVecCoeffs%abcof(i,l* (l+1)-l:l* (l+1)+l,1,natom,jsp),eigVecCoeffs%abcof(i,l* (l+1)-l:l* (l+1)+l,1,natom,jsp))
178 : ENDDO
179 879336 : regCharges%ener(l,n,jsp) = regCharges%ener(l,n,jsp) + (suma+sumb*usdus%ddn(l,n,jsp))/atoms%neq(n)*we(i)*eig(i)
180 947640 : regCharges%sqal(l,n,jsp) = regCharges%sqal(l,n,jsp) + (suma+sumb*usdus%ddn(l,n,jsp))/atoms%neq(n)*we(i)
181 : ENDDO
182 : ELSE
183 784 : n_dos=banddos%map_atomtype(n)
184 : !data present in dos-type
185 32656 : DO i = (skip_t+1),noccbd
186 31872 : regCharges%ener(l,n,jsp) = regCharges%ener(l,n,jsp) + dos%qal(l,n_dos,ev_list(i),ikpt,jsp)*we(i)*eig(i)
187 32656 : regCharges%sqal(l,n,jsp) = regCharges%sqal(l,n,jsp) + dos%qal(l,n_dos,ev_list(i),ikpt,jsp)*we(i)
188 : ENDDO
189 : ENDIF
190 : ENDDO
191 : ENDDO
192 :
193 : !---> initialize qlo
194 :
195 630374 : qlo(:,:,:,:) = 0.0
196 370138 : qaclo(:,:,:) = 0.0
197 370138 : qbclo(:,:,:) = 0.0
198 :
199 : !---> density for each local orbital and occupied state
200 :
201 25290 : DO ntyp=1,atoms%ntype
202 42866 : DO nn = 1,atoms%neq(ntyp)
203 17576 : natom = atoms%firstAtom(ntyp) - 1
204 17576 : natom = natom + nn
205 52428 : DO lo = 1,atoms%nlo(ntyp)
206 17580 : l = atoms%llo(lo,ntyp)
207 17580 : ll1 = l* (l+1)
208 59000 : DO m = -l,l
209 41420 : lm = ll1 + m
210 698783 : DO i = 1,ne
211 : qbclo(i,lo,ntyp) = qbclo(i,lo,ntyp) +REAL(&
212 : eigVecCoeffs%abcof(i,lm,1,natom,jsp)*CONJG(eigVecCoeffs%ccof(m,i,lo,natom,jsp))+&
213 639783 : eigVecCoeffs%ccof(m,i,lo,natom,jsp)*CONJG(eigVecCoeffs%abcof(i,lm,1,natom,jsp)) )
214 : qaclo(i,lo,ntyp) = qaclo(i,lo,ntyp) + REAL(&
215 : eigVecCoeffs%abcof(i,lm,0,natom,jsp)*CONJG(eigVecCoeffs%ccof(m,i,lo,natom,jsp))+&
216 681203 : eigVecCoeffs%ccof(m,i,lo,natom,jsp)*CONJG(eigVecCoeffs%abcof(i,lm,0,natom,jsp)) )
217 : ENDDO
218 : ENDDO
219 65336 : DO lop = 1,atoms%nlo(ntyp)
220 47760 : IF (atoms%llo(lop,ntyp).EQ.l) THEN
221 59160 : DO m = -l,l
222 704835 : DO i = 1,ne
223 : qlo(i,lop,lo,ntyp) = qlo(i,lop,lo,ntyp) + REAL(&
224 687215 : CONJG(eigVecCoeffs%ccof(m,i,lop,natom,jsp))*eigVecCoeffs%ccof(m,i,lo,natom,jsp))
225 : ENDDO
226 : ENDDO
227 : ENDIF
228 : ENDDO
229 : ENDDO
230 : ENDDO
231 : ENDDO
232 :
233 : !---> perform brillouin zone integration and sum over bands
234 :
235 25290 : DO ntyp = 1,atoms%ntype
236 17272 : n_dos=banddos%map_atomtype(ntyp)
237 42624 : DO lo = 1,atoms%nlo(ntyp)
238 17334 : l = atoms%llo(lo,ntyp)
239 : ! llo > 3 used for unoccupied states only
240 34606 : IF(l.LE.3) THEN
241 283474 : DO i = 1,ne
242 266140 : IF (n_dos>0) THEN
243 : dos%qal(l,n_dos,ev_list(i),ikpt,jsp)= dos%qal(l,n_dos,ev_list(i),ikpt,jsp) + ( 1.0/atoms%neq(ntyp) ) * &
244 5808 : (qaclo(i,lo,ntyp)*usdus%uulon(lo,ntyp,jsp)+qbclo(i,lo,ntyp)*usdus%dulon(lo,ntyp,jsp))
245 : ENDIF
246 : dos%qTot(ev_list(i),ikpt,jsp) = dos%qTot(ev_list(i),ikpt,jsp) + &
247 283474 : qaclo(i,lo,ntyp)*usdus%uulon(lo,ntyp,jsp)+qbclo(i,lo,ntyp)*usdus%dulon(lo,ntyp,jsp)
248 : END DO
249 47124 : DO lop = 1,atoms%nlo(ntyp)
250 47124 : IF (atoms%llo(lop,ntyp).EQ.l) THEN
251 285478 : DO i = 1,ne
252 268104 : regCharges%enerlo(lo,ntyp,jsp) = regCharges%enerlo(lo,ntyp,jsp) +qlo(i,lop,lo,ntyp)*we(i)*eig(i)
253 268104 : regCharges%sqlo(lo,ntyp,jsp) = regCharges%sqlo(lo,ntyp,jsp) + qlo(i,lop,lo,ntyp)*we(i)
254 268104 : IF (n_dos>0) THEN
255 : dos%qal(l,n_dos,ev_list(i),ikpt,jsp)= dos%qal(l,n_dos,ev_list(i),ikpt,jsp) + ( 1.0/atoms%neq(ntyp) ) *&
256 5808 : qlo(i,lop,lo,ntyp)*usdus%uloulopn(lop,lo,ntyp,jsp)
257 : END IF
258 285478 : dos%qTot(ev_list(i),ikpt,jsp) = dos%qTot(ev_list(i),ikpt,jsp) + qlo(i,lop,lo,ntyp)*usdus%uloulopn(lop,lo,ntyp,jsp)
259 : ENDDO
260 : ENDIF
261 : ENDDO
262 : ELSE
263 0 : DO i = 1,ne
264 : dos%qTot(ev_list(i),ikpt,jsp) = dos%qTot(ev_list(i),ikpt,jsp) + &
265 0 : qaclo(i,lo,ntyp)*usdus%uulon(lo,ntyp,jsp)+qbclo(i,lo,ntyp)*usdus%dulon(lo,ntyp,jsp)
266 : END DO
267 0 : DO lop = 1,atoms%nlo(ntyp)
268 0 : IF (atoms%llo(lop,ntyp).EQ.l) THEN
269 0 : DO i = 1,ne
270 0 : dos%qTot(ev_list(i),ikpt,jsp) = dos%qTot(ev_list(i),ikpt,jsp) + qlo(i,lop,lo,ntyp)*usdus%uloulopn(lop,lo,ntyp,jsp)
271 : ENDDO
272 : ENDIF
273 : ENDDO
274 : END IF
275 : END DO
276 : END DO
277 :
278 8018 : END SUBROUTINE eparas
279 : END MODULE m_eparas
|
