Line data Source code
1 : !--------------------------------------------------------------------------------
2 : ! Copyright (c) 2018 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_metagga
7 : PUBLIC :: calc_EnergyDen
8 : PRIVATE :: calc_EnergyDen_auxillary_weights, &
9 : calc_kinEnergyDen_pw, &
10 : calc_kinEnergyDen_mt
11 :
12 : type t_RS_potden
13 : REAL, ALLOCATABLE :: is(:,:), mt(:,:)
14 : end type t_RS_potden
15 :
16 : CONTAINS
17 0 : SUBROUTINE calc_kinEnergyDen_pw(EnergyDen_rs, vTot_rs, den_rs, kinEnergyDen_RS)
18 : USE m_juDFT_stop
19 : !use m_cdngen
20 : IMPLICIT NONE
21 : REAL, INTENT(in) :: den_RS(:,:), EnergyDen_RS(:,:), vTot_RS(:,:)
22 : REAL, INTENT(inout), allocatable :: kinEnergyDen_RS(:,:)
23 : #ifdef CPP_LIBXC
24 : REAL, PARAMETER :: eps = 1e-15
25 :
26 0 : kinEnergyDen_RS = EnergyDen_RS - vTot_RS * den_RS
27 : #else
28 : CALL juDFT_error("MetaGGA require LibXC",hint="compile Fleur with LibXC (e.g. by giving '-external libxc' to ./configure")
29 : #endif
30 0 : END SUBROUTINE calc_kinEnergyDen_pw
31 :
32 0 : SUBROUTINE calc_kinEnergyDen_mt(EnergyDen_RS, vTot_rs, vTot0_rs, core_den_rs, val_den_rs, &
33 0 : kinEnergyDen_RS)
34 : USE m_juDFT_stop
35 : USE m_juDFT_string
36 : implicit none
37 : REAL, INTENT(in) :: EnergyDen_RS(:,:), vTot_rs(:,:), vTot0_rs(:,:), core_den_rs(:,:), val_den_rs(:,:)
38 : REAL, INTENT(inout) :: kinEnergyDen_RS(:,:)
39 :
40 : #ifdef CPP_LIBXC
41 0 : kinEnergyDen_RS = EnergyDen_RS - (vTot0_rs * core_den_rs + vTot_rs * val_den_rs)
42 : #else
43 : CALL juDFT_error("MetaGGA require LibXC",hint="compile Fleur with LibXC (e.g. by giving '-external libxc' to ./configure")
44 : #endif
45 0 : END SUBROUTINE calc_kinEnergyDen_mt
46 :
47 :
48 0 : SUBROUTINE calc_EnergyDen(eig_id, mpi, kpts, noco, input, banddos, cell, atoms, enpara, stars, &
49 : vacuum, DIMENSION, sphhar, sym, vTot, oneD, results, EnergyDen)
50 : ! calculates the energy density
51 : ! EnergyDen = \sum_i n_i(r) \varepsilon_i
52 : ! where n_i(r) is the one-particle density
53 : ! and \varepsilon_i are the eigenenergies
54 :
55 : USE m_types_setup
56 : USE m_types_potden
57 : USE m_types_kpts
58 : USE m_types_mpi
59 : USE m_types_enpara
60 : USE m_types_misc
61 : USE m_types_regionCharges
62 : USE m_types_dos
63 : USE m_types_cdnval
64 : USE m_cdnval
65 :
66 : IMPLICIT NONE
67 :
68 : INTEGER, INTENT(in) :: eig_id
69 : TYPE(t_mpi), INTENT(in) :: mpi
70 : TYPE(t_kpts), INTENT(in) :: kpts
71 : TYPE(t_noco), INTENT(in) :: noco
72 : TYPE(t_input), INTENT(in) :: input
73 : TYPE(t_banddos), INTENT(in) :: banddos
74 : TYPE(t_cell), INTENT(in) :: cell
75 : TYPE(t_atoms), INTENT(in) :: atoms
76 : TYPE(t_enpara), INTENT(in) :: enpara
77 : TYPE(t_stars), INTENT(in) :: stars
78 : TYPE(t_vacuum), INTENT(in) :: vacuum
79 : TYPE(t_dimension), INTENT(in) :: DIMENSION
80 : TYPE(t_sphhar), INTENT(in) :: sphhar
81 : TYPE(t_sym), INTENT(in) :: sym
82 : TYPE(t_potden), INTENT(in) :: vTot
83 : TYPE(t_oneD), INTENT(in) :: oneD
84 : TYPE(t_results), INTENT(in) :: results
85 : TYPE(t_potden), INTENT(inout) :: EnergyDen
86 :
87 : ! local
88 : INTEGER :: jspin
89 :
90 0 : TYPE(t_regionCharges) :: regCharges
91 0 : TYPE(t_dos) :: dos
92 0 : TYPE(t_moments) :: moments
93 0 : TYPE(t_results) :: tmp_results
94 0 : TYPE(t_cdnvalJob) :: cdnvalJob
95 : TYPE(t_potden) :: aux_den, real_den
96 :
97 :
98 0 : CALL regCharges%init(input, atoms)
99 0 : CALL dos%init(input, atoms, DIMENSION, kpts, vacuum)
100 0 : CALL moments%init(input, atoms)
101 0 : tmp_results = results
102 :
103 0 : DO jspin = 1,input%jspins
104 0 : CALL cdnvalJob%init(mpi,input,kpts,noco,results,jspin)
105 :
106 :
107 : ! replace brillouin weights with auxillary weights
108 0 : CALL calc_EnergyDen_auxillary_weights(eig_id, kpts, jspin, cdnvalJob%weights)
109 :
110 : CALL cdnval(eig_id, mpi, kpts, jspin, noco, input, banddos, cell, atoms, &
111 : enpara, stars, vacuum, DIMENSION, sphhar, sym, vTot, oneD, cdnvalJob, &
112 0 : EnergyDen, regCharges, dos, tmp_results, moments)
113 : ENDDO
114 :
115 0 : END SUBROUTINE calc_EnergyDen
116 :
117 0 : SUBROUTINE calc_EnergyDen_auxillary_weights(eig_id, kpts, jspin, f_ik)
118 : USE m_types_kpts
119 : USE m_eig66_io
120 : IMPLICIT NONE
121 : ! calculates new (auxillary-)weights as
122 : ! f_iks = w_iks * E_iks
123 : !, where f_iks are the new (auxillary-)weights
124 : ! w_iks are the weights used in brillouin zone integration
125 : ! E_iks are the eigen energies
126 :
127 : INTEGER, INTENT(in) :: eig_id
128 : INTEGER, INTENT(in) :: jspin
129 : TYPE(t_kpts), INTENT(in) :: kpts
130 : REAL, INTENT(inout) :: f_ik(:,:) ! f_ik(band_idx, kpt_idx)
131 :
132 : ! local vars
133 0 : REAL :: e_i(SIZE(f_ik,dim=1))
134 : INTEGER :: ikpt
135 :
136 0 : DO ikpt = 1,kpts%nkpt
137 0 : CALL read_eig(eig_id,ikpt,jspin, eig=e_i)
138 0 : f_ik(:,ikpt) = f_ik(:,ikpt) * e_i
139 : ENDDO
140 0 : END SUBROUTINE calc_EnergyDen_auxillary_weights
141 :
142 0 : subroutine set_zPrime(dim_idx, zMat, kpt, lapw, cell, zPrime)
143 : USE m_types
144 : USE m_constants
145 : implicit none
146 : INTEGER, intent(in) :: dim_idx
147 : TYPE (t_mat), intent(in) :: zMat
148 : REAL, intent(in) :: kpt(3)
149 : TYPE(t_lapw), intent(in) :: lapw
150 : TYPE(t_cell), intent(in) :: cell
151 : TYPE (t_mat) :: zPrime
152 :
153 : REAL :: k_plus_g(3), fac
154 : INTEGER :: basis_idx
155 :
156 0 : call zPrime%free()
157 0 : call zPrime%init(zMat)
158 :
159 0 : do basis_idx = 1,size(lapw%gvec,dim=2)
160 0 : k_plus_g = kpt + lapw%gvec(:,basis_idx,1)
161 0 : k_plus_g = internal_to_rez(cell, k_plus_g)
162 :
163 0 : fac = k_plus_g(dim_idx)
164 0 : if(zPrime%l_real) then
165 0 : zPrime%data_r(basis_idx,:) = fac * zMat%data_r(basis_idx,:)
166 : else
167 0 : zPrime%data_c(basis_idx,:) = ImagUnit * fac * zMat%data_c(basis_idx,:)
168 : endif
169 : enddo
170 0 : end subroutine set_zPrime
171 :
172 0 : function internal_to_rez(cell, vec) result(res)
173 : use m_types
174 : implicit none
175 : type(t_cell), intent(in) :: cell
176 : real, intent(in) :: vec(3)
177 : real :: res(3)
178 :
179 0 : res = matmul(transpose(cell%bmat), vec)
180 0 : end function internal_to_rez
181 :
182 0 : subroutine undo_vgen_finalize(vtot, atoms, noco, stars)
183 : use m_types
184 : use m_constants
185 : use m_judft
186 : implicit none
187 : TYPE(t_potden), intent(inout) :: vtot
188 : type(t_atoms), intent(in) :: atoms
189 : type(t_noco), intent(in) :: noco
190 : type(t_stars), intent(in) :: stars
191 :
192 : integer :: js, n, st
193 :
194 0 : do js = 1,size(vtot%mt,4)
195 0 : do n = 1,atoms%ntype
196 : vTot%mt(:atoms%jri(n),0,n,js) = vtot%mt(:atoms%jri(n),0,n,js) &
197 0 : / (atoms%rmsh(:atoms%jri(n),n) / sfp_const )
198 : enddo
199 : enddo
200 :
201 0 : if(.not. noco%l_noco) then
202 0 : do js=1,size(vtot%pw_w,2)
203 0 : do st=1,stars%ng3
204 0 : vTot%pw_w(st,js) = vTot%pw_w(st,js) * stars%nstr(st)
205 : enddo
206 : enddo
207 : else
208 0 : call juDFT_error("undo vgen_finalize not implemented for noco")
209 : endif
210 0 : end subroutine undo_vgen_finalize
211 :
212 0 : subroutine set_kinED(mpi, sphhar, atoms, sym, core_den, val_den, xcpot, &
213 : input, noco, stars, cell, den, EnergyDen, vTot)
214 : use m_types
215 : implicit none
216 : TYPE(t_mpi),INTENT(IN) :: mpi
217 : TYPE(t_sphhar),INTENT(IN) :: sphhar
218 : TYPE(t_atoms),INTENT(IN) :: atoms
219 : TYPE(t_sym), INTENT(IN) :: sym
220 : TYPE(t_potden),INTENT(IN) :: core_den, val_den
221 : CLASS(t_xcpot),INTENT(INOUT) :: xcpot
222 : TYPE(t_input),INTENT(IN) :: input
223 : TYPE(t_noco),INTENT(IN) :: noco
224 : TYPE(t_stars),INTENT(IN) :: stars
225 : TYPE(t_cell),INTENT(IN) :: cell
226 : TYPE(t_potden),INTENT(IN) :: den, EnergyDen, vTot
227 :
228 0 : TYPE(t_potden) :: vTot_corrected
229 :
230 0 : call vTot_corrected%copyPotDen(vTot)
231 0 : call undo_vgen_finalize(vTot_corrected, atoms, noco, stars)
232 :
233 0 : call set_kinED_is(xcpot, input, noco, stars, sym, cell, den, EnergyDen, vTot_corrected)
234 : call set_kinED_mt(mpi, sphhar, atoms, sym, core_den, val_den, &
235 0 : xcpot, EnergyDen, input, vTot_corrected)
236 0 : end subroutine set_kinED
237 :
238 0 : subroutine set_kinED_is(xcpot, input, noco, stars, sym, cell, den, EnergyDen, vTot)
239 : use m_types
240 : use m_pw_tofrom_grid
241 : implicit none
242 : CLASS(t_xcpot),INTENT(INOUT) :: xcpot
243 : TYPE(t_input),INTENT(IN) :: input
244 : TYPE(t_noco),INTENT(IN) :: noco
245 : TYPE(t_stars),INTENT(IN) :: stars
246 : TYPE(t_sym), INTENT(IN) :: sym
247 : TYPE(t_cell),INTENT(IN) :: cell
248 : TYPE(t_potden),INTENT(IN) :: den, EnergyDen, vTot
249 :
250 : !local arrays
251 0 : REAL, ALLOCATABLE :: den_rs(:,:), ED_rs(:,:), vTot_rs(:,:)
252 0 : TYPE(t_gradients) :: tmp_grad
253 :
254 0 : CALL init_pw_grid(xcpot,stars,sym,cell)
255 :
256 : CALL pw_to_grid(xcpot, input%jspins, noco%l_noco, stars, &
257 : cell, EnergyDen%pw, tmp_grad, ED_rs)
258 : CALL pw_to_grid(xcpot, input%jspins, noco%l_noco, stars, &
259 : cell, vTot%pw, tmp_grad, vTot_rs)
260 : CALL pw_to_grid(xcpot, input%jspins, noco%l_noco, stars, &
261 : cell, den%pw, tmp_grad, den_rs)
262 :
263 0 : CALL finish_pw_grid()
264 :
265 0 : call calc_kinEnergyDen_pw(ED_rs, vTot_rs, den_rs, xcpot%kinED%is)
266 : !xcpot%kinED%is = ED_RS - vTot_RS * den_RS
267 0 : xcpot%kinED%set = .True.
268 0 : end subroutine set_kinED_is
269 :
270 0 : subroutine set_kinED_mt(mpi, sphhar, atoms, sym, core_den, val_den, &
271 : xcpot, EnergyDen, input, vTot)
272 : use m_types
273 : use m_mt_tofrom_grid
274 : implicit none
275 : TYPE(t_mpi),INTENT(IN) :: mpi
276 : TYPE(t_sphhar),INTENT(IN) :: sphhar
277 : TYPE(t_atoms),INTENT(IN) :: atoms
278 : TYPE(t_sym), INTENT(IN) :: sym
279 : TYPE(t_potden),INTENT(IN) :: core_den, val_den, EnergyDen, vTot
280 : CLASS(t_xcpot),INTENT(INOUT) :: xcpot
281 : TYPE(t_input),INTENT(IN) :: input
282 :
283 : INTEGER :: jr, loc_n, n, n_start, n_stride, cnt
284 0 : REAL,ALLOCATABLE :: vTot_mt(:,:,:), ED_rs(:,:), vTot_rs(:,:), vTot0_rs(:,:),&
285 0 : core_den_rs(:,:), val_den_rs(:,:)
286 0 : TYPE(t_gradients) :: tmp_grad
287 0 : TYPE(t_sphhar) :: tmp_sphhar
288 :
289 : #ifdef CPP_MPI
290 0 : n_start=mpi%irank+1
291 0 : n_stride=mpi%isize
292 : #else
293 : n_start=1
294 : n_stride=1
295 : #endif
296 0 : CALL init_mt_grid(input%jspins,atoms,sphhar,xcpot,sym)
297 0 : loc_n = 0
298 0 : allocate(ED_rs(atoms%nsp()*atoms%jmtd, input%jspins))
299 0 : allocate(vTot_rs, mold=ED_rs)
300 0 : allocate(vTot0_rs, mold=ED_rs)
301 0 : allocate(core_den_rs, mold=ED_rs)
302 0 : allocate(val_den_rs, mold=ED_rs)
303 :
304 : call xcpot%kinED%alloc_mt(atoms%nsp()*atoms%jmtd, input%jspins, &
305 0 : n_start, atoms%ntype, n_stride)
306 0 : loc_n = 0
307 0 : do n = n_start,atoms%ntype,n_stride
308 0 : loc_n = loc_n + 1
309 :
310 0 : if(.not. allocated(vTot_mt)) then
311 : allocate(vTot_mt(lbound(vTot%mt, dim=1):ubound(vTot%mt, dim=1),&
312 : lbound(vTot%mt, dim=2):ubound(vTot%mt, dim=2),&
313 0 : lbound(vTot%mt, dim=4):ubound(vTot%mt, dim=4)))
314 : endif
315 :
316 0 : do jr=1,atoms%jri(n)
317 0 : vTot_mt(jr,0:,:) = vTot%mt(jr,0:,n,:) * atoms%rmsh(jr,n)**2
318 : enddo
319 : CALL mt_to_grid(xcpot, input%jspins, atoms, sphhar, EnergyDen%mt(:, 0:, n, :), &
320 0 : n, tmp_grad, ED_rs)
321 : CALL mt_to_grid(xcpot, input%jspins, atoms, sphhar, vTot_mt(:,0:,:), &
322 0 : n, tmp_grad, vTot_rs)
323 :
324 0 : tmp_sphhar%nlhd = sphhar%nlhd
325 0 : tmp_sphhar%nlh = [(0, cnt=1,size(sphhar%nlh))]
326 :
327 : CALL mt_to_grid(xcpot, input%jspins, atoms, tmp_sphhar, vTot_mt(:,0:0,:), &
328 0 : n, tmp_grad, vTot0_rs)
329 : CALL mt_to_grid(xcpot, input%jspins, atoms, sphhar, &
330 0 : core_den%mt(:,0:,n,:), n, tmp_grad, core_den_rs)
331 : CALL mt_to_grid(xcpot, input%jspins, atoms, sphhar, &
332 0 : val_den%mt(:,0:,n,:), n, tmp_grad, val_den_rs)
333 :
334 : call calc_kinEnergyDen_mt(ED_RS, vTot_rs, vTot0_rs, core_den_rs, val_den_rs, &
335 0 : xcpot%kinED%mt(:,:,loc_n))
336 : !xcpot%kinED%mt(:,:,loc_n) = ED_RS - (vTot0_rs * core_den_rs + vTot_rs * val_den_rs)
337 : enddo
338 0 : xcpot%kinED%set = .True.
339 0 : CALL finish_mt_grid()
340 0 : end subroutine set_kinED_mt
341 : END MODULE m_metagga
|
