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 : MODULE m_doubleCounting
7 :
8 : USE m_constants
9 : USE m_types
10 : USE m_coulombPotential
11 :
12 : IMPLICIT NONE
13 :
14 : CONTAINS
15 :
16 160 : FUNCTION doubleCountingPot(density, ldau, U,J, l_spinoffd, l_mix,l_spinAvg, alpha, l_write) RESULT(Vdc)
17 :
18 : !------------------------------------------------------
19 : ! Calculate the Double Counting Correction in either
20 : ! the FLL or AMF limit
21 : ! If l_spinAvg is True the Double counting will be
22 : ! averaged over the spins
23 : !------------------------------------------------------
24 :
25 : COMPLEX, INTENT(IN) :: density(-lmaxU_const:,-lmaxU_const:,:)
26 : TYPE(t_utype), INTENT(IN) :: ldau !LDA+U information
27 : REAL, INTENT(IN) :: U, J
28 : LOGICAL, INTENT(IN) :: l_spinoffd
29 : LOGICAL, INTENT(IN) :: l_mix !Mix between FLL and AMF
30 : LOGICAL, INTENT(IN) :: l_spinAvg !Do we want a spin averaged double counting
31 : REAL, INTENT(IN) :: alpha
32 : LOGICAL, OPTIONAL, INTENT(IN) :: l_write
33 :
34 : COMPLEX :: Vdc(-lmaxU_const:lmaxU_const,-lmaxU_const:lmaxU_const,SIZE(density,3))
35 :
36 160 : COMPLEX, ALLOCATABLE :: modified_density(:,:,:)
37 : REAL :: charge, mag(3),mag_m(0:3), tmp, D, Vdcup,Vdcdn
38 : COMPLEX :: sigma(2,2,3), r21
39 : INTEGER :: spin_dim, ispin,m, spin1,spin2,mp
40 640 : type(t_nococonv) :: nococonv !Used only for the procedures on it
41 :
42 160 : sigma = cmplx_0
43 160 : sigma(1,2,1)=CMPLX(1.0,0.0)
44 160 : sigma(2,1,1)=CMPLX(1.0,0.0)
45 160 : sigma(1,2,2)=CMPLX(0.0,-1.0)
46 160 : sigma(2,1,2)=CMPLX(0.0,1.0)
47 160 : sigma(1,1,3)=CMPLX(1.0,0.0)
48 160 : sigma(2,2,3)=CMPLX(-1.0,0.0)
49 :
50 160 : spin_dim = SIZE(density,3)
51 160 : IF(.NOT.l_spinoffd) spin_dim = MIN(2,spin_dim)
52 :
53 160 : charge = 0.0
54 160 : mag = 0.0
55 :
56 880 : DO m = -ldau%l, ldau%l
57 720 : if (spin_dim==3) then
58 0 : r21 = density(m,m,3)
59 : else
60 720 : r21 = 0.0
61 : endif
62 : mag_m = nococonv%denmat_to_mag(real(density(m,m,1)),&
63 : real(density(m,m,min(2,spin_dim))),&
64 720 : r21)
65 :
66 720 : charge = charge + mag_m(0)
67 3040 : mag = mag + mag_m(1:)
68 : ENDDO
69 160 : IF(spin_dim == 1) then
70 104 : charge=charge/2.0
71 104 : mag = 0.0
72 : ENDIF
73 :
74 12472 : Vdc = cmplx_0
75 160 : IF(ldau%l_amf) THEN
76 0 : ALLOCATE(modified_density(-lmaxU_const:lmaxU_const, -lmaxU_const:lmaxU_const, SIZE(density,3)), source=cmplx_0)
77 0 : modified_density = cmplx_0
78 0 : D = real(2*(2*ldau%l+1))
79 :
80 0 : DO ispin = 1, spin_dim
81 :
82 0 : IF(ispin==3) THEN
83 : spin1 = 2
84 : spin2 = 1
85 : ELSE
86 0 : spin1 = ispin
87 0 : spin2 = ispin
88 : ENDIF
89 :
90 0 : DO m = -ldau%l, ldau%l
91 0 : IF(spin1==spin2) modified_density(m,m,ispin) = charge/D
92 0 : modified_density(m,m,ispin) = modified_density(m,m,ispin) + dot_product(mag,sigma(spin1,spin2,:))/D
93 : ENDDO
94 :
95 : ENDDO
96 0 : call coulombPotential(modified_density,ldau, MIN(2,SIZE(density,3)), l_spinoffd,Vdc,tmp)
97 : ELSE
98 376 : DO ispin = 1, spin_dim
99 :
100 216 : IF(ispin==3) THEN
101 : spin1 = 2
102 : spin2 = 1
103 : ELSE
104 216 : spin1 = ispin
105 216 : spin2 = ispin
106 : ENDIF
107 :
108 1400 : DO m = -ldau%l, ldau%l
109 1024 : IF(spin1 == spin2) Vdc(m,m,ispin) = Vdc(m,m,ispin) + U*(charge-0.5) - J*(charge/2.0-0.5)
110 4312 : Vdc(m,m,ispin) = Vdc(m,m,ispin) - J/2.0*dot_product(mag,sigma(spin1,spin2,:))
111 : ENDDO
112 : ENDDO
113 : ENDIF
114 :
115 :
116 160 : IF(PRESENT(l_write)) THEN
117 0 : IF(l_write) THEN
118 0 : WRITE(oUnit,"(/,A)") 'Double counting chemical potential:'
119 0 : IF(ldau%l_amf) THEN
120 0 : WRITE(oUnit,9040) 'AMF: ','spin-up','spin-dn','(up+dn)/2','up-dn'
121 : ELSE
122 0 : WRITE(oUnit,9040) 'FLL: ','spin-up','spin-dn','(up+dn)/2','up-dn'
123 : ENDIF
124 : 9040 FORMAT(TR3,A4,TR1,A7,TR3,A7,TR3,A9,TR3,A5)
125 0 : Vdcup = 0.0
126 0 : Vdcdn = 0.0
127 0 : do m = -ldau%l, ldau%l
128 0 : Vdcup = Vdcup + Vdc(m,m,1)/real(2*ldau%l+1)
129 0 : Vdcdn = Vdcdn + Vdc(m,m,min(2,spin_dim))/real(2*ldau%l+1)
130 : enddo
131 0 : WRITE(oUnit,9050) Vdcup,Vdcdn,(Vdcup+Vdcdn)/2.0,Vdcup-Vdcdn
132 : 9050 FORMAT(TR7,f8.4,TR2,f8.4,TR2,f8.4,TR4,f8.4)
133 : ENDIF
134 : ENDIF
135 :
136 160 : IF(l_spinAvg) THEN
137 0 : Vdc(:,:,1) = (Vdc(:,:,1)+Vdc(:,:,min(2,spin_dim)))/2.0
138 0 : if(spin_dim>1) Vdc(:,:,2) = Vdc(:,:,1)
139 0 : if(spin_dim==3) Vdc(:,:,3) = cmplx_0 !Is this right?
140 : ENDIF
141 :
142 160 : END FUNCTION doubleCountingPot
143 :
144 :
145 160 : REAL FUNCTION doubleCountingEnergy(density, ldau, U,J, l_spinoffd, l_mix,l_spinAvg, alpha, l_write)
146 :
147 : !------------------------------------------------------------
148 : ! Calculate the Double Counting Correction Energy in either
149 : ! the FLL or AMF limit
150 : !------------------------------------------------------------
151 :
152 : COMPLEX, INTENT(IN) :: density(-lmaxU_const:,-lmaxU_const:,:)
153 : TYPE(t_utype), INTENT(IN) :: ldau !LDA+U information
154 : REAL, INTENT(IN) :: U, J
155 : LOGICAL, INTENT(IN) :: l_spinoffd
156 : LOGICAL, INTENT(IN) :: l_mix !Mix between FLL and AMF
157 : LOGICAL, INTENT(IN) :: l_spinAvg !Do we want a spin averaged double counting
158 : REAL, INTENT(IN) :: alpha
159 : LOGICAL, OPTIONAL, INTENT(IN) :: l_write
160 :
161 : REAL :: charge, mag(3), mag_m(0:3), D
162 : INTEGER :: spin_dim, ispin,m, spin1,spin2
163 160 : COMPLEX :: tmp(-lmaxU_const:lmaxU_const,-lmaxU_const:lmaxU_const,SIZE(density,3)), r21
164 160 : COMPLEX, ALLOCATABLE :: modified_density(:,:,:)
165 : COMPLEX :: sigma(2,2,3)
166 640 : type(t_nococonv) :: nococonv !Used only for the procedures on it
167 :
168 160 : sigma = cmplx_0
169 160 : sigma(1,2,1)=CMPLX(1.0,0.0)
170 160 : sigma(2,1,1)=CMPLX(1.0,0.0)
171 160 : sigma(1,2,2)=CMPLX(0.0,-1.0)
172 160 : sigma(2,1,2)=CMPLX(0.0,1.0)
173 160 : sigma(1,1,3)=CMPLX(1.0,0.0)
174 160 : sigma(2,2,3)=CMPLX(-1.0,0.0)
175 :
176 160 : spin_dim = SIZE(density,3)
177 160 : IF(.NOT.l_spinoffd) spin_dim = MIN(2,spin_dim)
178 :
179 160 : charge = 0.0
180 160 : mag = 0.0
181 880 : DO m = -ldau%l, ldau%l
182 720 : if (spin_dim==3) then
183 0 : r21 = density(m,m,3)
184 : else
185 720 : r21 = 0.0
186 : endif
187 : mag_m = nococonv%denmat_to_mag(real(density(m,m,1)),&
188 : real(density(m,m,min(2,spin_dim))),&
189 720 : r21)
190 :
191 720 : charge = charge + mag_m(0)
192 3040 : mag = mag + mag_m(1:)
193 : ENDDO
194 :
195 160 : IF(spin_dim == 1) then
196 104 : charge=charge/2.0
197 104 : mag = 0.
198 : endif
199 160 : if (l_spinAvg) mag = 0.
200 :
201 160 : doubleCountingEnergy = 0.0
202 160 : IF(ldau%l_amf) THEN
203 0 : ALLOCATE(modified_density(-lmaxU_const:lmaxU_const, -lmaxU_const:lmaxU_const, SIZE(density,3)), source=cmplx_0)
204 0 : modified_density = cmplx_0
205 0 : D = real(2*(2*ldau%l+1))
206 :
207 0 : DO ispin = 1, spin_dim
208 :
209 0 : IF(ispin==3) THEN
210 : spin1 = 2
211 : spin2 = 1
212 : ELSE
213 0 : spin1 = ispin
214 0 : spin2 = ispin
215 : ENDIF
216 :
217 0 : DO m = -ldau%l, ldau%l
218 0 : IF(spin1==spin2) modified_density(m,m,ispin) = charge/D
219 0 : modified_density(m,m,ispin) = modified_density(m,m,ispin) + dot_product(mag,sigma(spin1,spin2,:))/D
220 : ENDDO
221 :
222 : ENDDO
223 0 : call coulombPotential(density-modified_density,ldau, MIN(2,SIZE(density,3)), l_spinoffd,tmp,doubleCountingEnergy)
224 : ELSE
225 640 : doubleCountingEnergy = U/2*charge*(charge-1) -J/2*charge*(charge/2-1)-J*dot_product(mag,mag)/4
226 : ENDIF
227 :
228 160 : END FUNCTION doubleCountingEnergy
229 :
230 0 : FUNCTION doubleCountingMixFactor(mmpmat, l, charge) Result(alpha)
231 :
232 : !---------------------------------------------
233 : ! Calculate the mixing factor between FLL/AMF
234 : !---------------------------------------------
235 :
236 : COMPLEX, INTENT(IN) :: mmpmat(-lmaxU_const:, -lmaxU_const:, :)
237 : INTEGER, INTENT(IN) :: l
238 : REAL, INTENT(IN) :: charge
239 :
240 : REAL :: alpha
241 :
242 0 : alpha = 0.0
243 :
244 0 : END FUNCTION
245 :
246 : END MODULE m_doubleCounting
|
