(git:374b731)
Loading...
Searching...
No Matches
atom_xc.F
Go to the documentation of this file.
1!--------------------------------------------------------------------------------------------------!
2! CP2K: A general program to perform molecular dynamics simulations !
3! Copyright 2000-2024 CP2K developers group <https://cp2k.org> !
4! !
5! SPDX-License-Identifier: GPL-2.0-or-later !
6!--------------------------------------------------------------------------------------------------!
7
8! **************************************************************************************************
9!> \brief routines that build the integrals of the Vxc potential calculated
10!> for the atomic code
11! **************************************************************************************************
12MODULE atom_xc
13
14 USE atom_types, ONLY: atom_basis_type,&
15 atom_type,&
16 gth_pseudo,&
17 opmat_type,&
18 sgp_pseudo,&
19 upf_pseudo
20 USE atom_utils, ONLY: atom_core_density,&
21 atom_density,&
22 coulomb_potential_numeric,&
23 integrate_grid,&
24 numpot_matrix
25 USE cp_files, ONLY: close_file,&
26 get_unit_number,&
27 open_file
28 USE input_constants, ONLY: xc_none
29 USE input_section_types, ONLY: section_vals_get_subs_vals,&
30 section_vals_type,&
31 section_vals_val_get
32 USE kinds, ONLY: dp
33 USE mathconstants, ONLY: fourpi,&
34 pi
35 USE xc_atom, ONLY: xc_rho_set_atom_update
36 USE xc_derivative_desc, ONLY: &
37 deriv_norm_drho, deriv_norm_drhoa, deriv_norm_drhob, deriv_rho, deriv_rhoa, deriv_rhob, &
38 deriv_tau, deriv_tau_a, deriv_tau_b
39 USE xc_derivative_set_types, ONLY: xc_derivative_set_type,&
40 xc_dset_create,&
41 xc_dset_get_derivative,&
42 xc_dset_release,&
43 xc_dset_zero_all
44 USE xc_derivative_types, ONLY: xc_derivative_get,&
45 xc_derivative_type
46 USE xc_derivatives, ONLY: xc_functionals_eval,&
47 xc_functionals_get_needs
48 USE xc_rho_cflags_types, ONLY: xc_rho_cflags_type
49 USE xc_rho_set_types, ONLY: xc_rho_set_create,&
50 xc_rho_set_release,&
51 xc_rho_set_type
52#include "./base/base_uses.f90"
53
54 IMPLICIT NONE
55
56 PRIVATE
57
58 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'atom_xc'
59
60! ZMP added coulomb_potential_numeric
61! ZMP added public subroutines calculate_atom_zmp, calculate_atom_ext_vxc
62 PUBLIC :: calculate_atom_vxc_lda, calculate_atom_vxc_lsd, &
63 calculate_atom_zmp, calculate_atom_ext_vxc
64 PUBLIC :: atom_vxc_lda, atom_vxc_lsd, atom_dpot_lda
65
66CONTAINS
67
68! **************************************************************************************************
69!> \brief ZMP subroutine for the calculation of the effective constraint
70!> potential in the atomic code.
71!> \param ext_density external electron density
72!> \param atom information about the atomic kind
73!> \param lprint save exchange-correlation potential to the file 'linear_potentials.dat'
74!> \param xcmat exchange-correlation potential matrix
75!> \author D. Varsano [daniele.varsano@nano.cnr.it]
76! **************************************************************************************************
77 SUBROUTINE calculate_atom_zmp(ext_density, atom, lprint, xcmat)
78 REAL(kind=dp), DIMENSION(:), INTENT(INOUT) :: ext_density
79 TYPE(atom_type), INTENT(INOUT) :: atom
80 LOGICAL :: lprint
81 TYPE(opmat_type), OPTIONAL, POINTER :: xcmat
82
83 CHARACTER(LEN=*), PARAMETER :: routinen = 'calculate_atom_zmp'
84
85 INTEGER :: extunit, handle, ir, nr, z
86 REAL(kind=dp) :: int1, int2
87 REAL(kind=dp), DIMENSION(:), POINTER :: deltarho, rho_dum, vxc, vxc1, vxc2
88 REAL(kind=dp), DIMENSION(:, :), POINTER :: rho
89
90 CALL timeset(routinen, handle)
91
92 nr = atom%basis%grid%nr
93 z = atom%z
94 ALLOCATE (rho(nr, 1), vxc(nr), vxc1(nr), vxc2(nr), rho_dum(nr), deltarho(nr))
95 CALL atom_density(rho(:, 1), atom%orbitals%pmat, atom%basis, atom%state%maxl_occ, typ="RHO")
96 !Vxc1
97 int1 = integrate_grid(ext_density, atom%basis%grid)
98 int1 = fourpi*int1
99 CALL coulomb_potential_numeric(vxc1, rho(:, 1), atom%basis%grid)
100
101 vxc1 = -vxc1/z
102
103 !Vxc2
104 rho_dum = rho(:, 1)*int1/z
105 deltarho = rho_dum - ext_density
106 int2 = integrate_grid(deltarho, atom%basis%grid)
107 CALL coulomb_potential_numeric(vxc2, deltarho, atom%basis%grid)
108
109 vxc2 = vxc2*atom%lambda
110
111 !Vxc
112 vxc = vxc1 + vxc2
113
114 atom%energy%exc = fourpi*integrate_grid(vxc, rho(:, 1), atom%basis%grid)
115 atom%rho_diff_integral = fourpi*int2
116
117 IF (lprint) THEN
118 extunit = get_unit_number()
119 CALL open_file(file_name="linear_potentials.dat", file_status="UNKNOWN", &
120 file_form="FORMATTED", file_action="WRITE", &
121 unit_number=extunit)
122
123 WRITE (extunit, fmt='("#",T10,"R[bohr]",T36,"Rho[au]",T61,"DRho[au]",T86,"V_XC[au]",T111,"V_XC1[au]",T136,"V_XC2[au]")')
124 DO ir = 1, nr
125 WRITE (extunit, fmt='(T1,E24.15,T26,E24.15,T51,E24.15,T76,E24.15,T101,E24.15,T126,E24.15)') &
126 atom%basis%grid%rad(ir), rho(ir, 1), deltarho(ir), vxc(ir), vxc1(ir), vxc2(ir)
127 END DO
128 CALL close_file(unit_number=extunit)
129 END IF
130
131 IF (PRESENT(xcmat)) CALL numpot_matrix(xcmat%op, vxc, atom%basis, 0)
132
133 DEALLOCATE (rho, vxc, vxc1, vxc2, rho_dum, deltarho)
134
135 CALL timestop(handle)
136
137 END SUBROUTINE calculate_atom_zmp
138
139! **************************************************************************************************
140!> \brief ZMP subroutine for the scf external density from the external v_xc.
141!> \param vxc exchange-correlation potential
142!> \param atom information about the atomic kind
143!> \param lprint save exchange-correlation potential to the file 'linear_potentials.dat'
144!> \param xcmat exchange-correlation potential matrix
145!> \author D. Varsano [daniele.varsano@nano.cnr.it]
146! **************************************************************************************************
147 SUBROUTINE calculate_atom_ext_vxc(vxc, atom, lprint, xcmat)
148 REAL(kind=dp), DIMENSION(:), INTENT(IN) :: vxc
149 TYPE(atom_type), INTENT(INOUT) :: atom
150 LOGICAL, INTENT(in) :: lprint
151 TYPE(opmat_type), INTENT(inout), OPTIONAL, POINTER :: xcmat
152
153 CHARACTER(LEN=*), PARAMETER :: routinen = 'calculate_atom_ext_vxc'
154
155 INTEGER :: extunit, handle, ir, nr
156 REAL(kind=dp), DIMENSION(:, :), POINTER :: rho
157
158 CALL timeset(routinen, handle)
159 nr = atom%basis%grid%nr
160
161 ALLOCATE (rho(nr, 1))
162
163 CALL atom_density(rho(:, 1), atom%orbitals%pmat, atom%basis, atom%state%maxl_occ, typ="RHO")
164
165 IF (lprint) THEN
166 extunit = get_unit_number()
167 CALL open_file(file_name="linear_potentials.dat", file_status="UNKNOWN", &
168 file_form="FORMATTED", file_action="WRITE", &
169 unit_number=extunit)
170
171 WRITE (extunit, fmt='("#",T10,"R[bohr]",T36,"Rho[au]",T61,"V_XC[au]")')
172 DO ir = 1, nr
173 WRITE (extunit, fmt='(T1,E24.15,T26,E24.15,T51,E24.15)') &
174 atom%basis%grid%rad(ir), rho(ir, 1), vxc(ir)
175 END DO
176 CALL close_file(unit_number=extunit)
177 END IF
178
179 atom%energy%exc = fourpi*integrate_grid(vxc, rho(:, 1), atom%basis%grid)
180 IF (PRESENT(xcmat)) CALL numpot_matrix(xcmat%op, vxc, atom%basis, 0)
181
182 DEALLOCATE (rho)
183
184 CALL timestop(handle)
185
186 END SUBROUTINE calculate_atom_ext_vxc
187
188! **************************************************************************************************
189!> \brief Calculate atomic exchange-correlation potential in spin-restricted case.
190!> \param xcmat exchange-correlation potential expressed in the atomic basis set
191!> \param atom information about the atomic kind
192!> \param xc_section XC input section
193!> \par History
194!> * 03.2010 extension of GTH pseudo-potential definition [Juerg Hutter]
195!> * 02.2010 renamed to calculate_atom_vxc_lda() [Juerg Hutter]
196!> * 08.2008 created [Juerg Hutter]
197! **************************************************************************************************
198 SUBROUTINE calculate_atom_vxc_lda(xcmat, atom, xc_section)
199 TYPE(opmat_type), POINTER :: xcmat
200 TYPE(atom_type), INTENT(INOUT) :: atom
201 TYPE(section_vals_type), POINTER :: xc_section
202
203 CHARACTER(LEN=*), PARAMETER :: routinen = 'calculate_atom_vxc_lda'
204
205 INTEGER :: deriv_order, handle, i, l, myfun, n1, &
206 n2, n3, nr, nspins, unit_nr
207 INTEGER, DIMENSION(2, 3) :: bounds
208 LOGICAL :: lsd, nlcc
209 REAL(kind=dp) :: density_cut, gradient_cut, tau_cut
210 REAL(kind=dp), DIMENSION(:), POINTER :: exc, vxc
211 REAL(kind=dp), DIMENSION(:, :), POINTER :: drho, lap, rho, tau
212 REAL(kind=dp), DIMENSION(:, :, :), POINTER :: taumat, xcpot
213 TYPE(section_vals_type), POINTER :: xc_fun_section
214 TYPE(xc_derivative_set_type) :: deriv_set
215 TYPE(xc_derivative_type), POINTER :: deriv
216 TYPE(xc_rho_cflags_type) :: needs
217 TYPE(xc_rho_set_type) :: rho_set
218
219 CALL timeset(routinen, handle)
220
221 nlcc = .false.
222 IF (atom%potential%ppot_type == gth_pseudo) THEN
223 nlcc = atom%potential%gth_pot%nlcc
224 ELSE IF (atom%potential%ppot_type == upf_pseudo) THEN
225 nlcc = atom%potential%upf_pot%core_correction
226 ELSE IF (atom%potential%ppot_type == sgp_pseudo) THEN
227 nlcc = atom%potential%sgp_pot%has_nlcc
228 END IF
229
230 IF (ASSOCIATED(xc_section)) THEN
231 xc_fun_section => section_vals_get_subs_vals(xc_section, "XC_FUNCTIONAL")
232 CALL section_vals_val_get(xc_fun_section, "_SECTION_PARAMETERS_", i_val=myfun)
233
234 IF (myfun == xc_none) THEN
235 atom%energy%exc = 0._dp
236 ELSE
237 CALL section_vals_val_get(xc_section, "DENSITY_CUTOFF", r_val=density_cut)
238 CALL section_vals_val_get(xc_section, "GRADIENT_CUTOFF", r_val=gradient_cut)
239 CALL section_vals_val_get(xc_section, "TAU_CUTOFF", r_val=tau_cut)
240
241 lsd = .false.
242 nspins = 1
243 needs = xc_functionals_get_needs(xc_fun_section, lsd=lsd, calc_potential=.false.)
244
245 ! Prepare the structures needed to calculate and store the xc derivatives
246
247 ! Array dimension: here anly one dimensional arrays are used,
248 ! i.e. only the first column of deriv_data is read.
249 ! The other to dimensions are set to size equal 1
250 nr = atom%basis%grid%nr
251 bounds(1:2, 1:3) = 1
252 bounds(2, 1) = nr
253
254 ! create a place where to put the derivatives
255 CALL xc_dset_create(deriv_set, local_bounds=bounds)
256 ! create the place where to store the argument for the functionals
257 CALL xc_rho_set_create(rho_set, bounds, rho_cutoff=density_cut, &
258 drho_cutoff=gradient_cut, tau_cutoff=tau_cut)
259 ! allocate the required 3d arrays where to store rho and drho
260 CALL xc_rho_set_atom_update(rho_set, needs, nspins, bounds)
261
262 NULLIFY (rho, drho, tau)
263 IF (needs%rho) THEN
264 ALLOCATE (rho(nr, 1))
265 CALL atom_density(rho(:, 1), atom%orbitals%pmat, atom%basis, atom%state%maxl_occ, typ="RHO")
266 IF (nlcc) THEN
267 CALL atom_core_density(rho(:, 1), atom%potential, typ="RHO", rr=atom%basis%grid%rad)
268 END IF
269 END IF
270 IF (needs%norm_drho) THEN
271 ALLOCATE (drho(nr, 1))
272 CALL atom_density(drho(:, 1), atom%orbitals%pmat, atom%basis, atom%state%maxl_occ, typ="DER")
273 IF (nlcc) THEN
274 CALL atom_core_density(drho(:, 1), atom%potential, typ="DER", rr=atom%basis%grid%rad)
275 END IF
276 END IF
277 IF (needs%tau) THEN
278 ALLOCATE (tau(nr, 1))
279 CALL atom_density(tau(:, 1), atom%orbitals%pmat, atom%basis, atom%state%maxl_occ, &
280 typ="KIN", rr=atom%basis%grid%rad2)
281 END IF
282 IF (needs%laplace_rho) THEN
283 ALLOCATE (lap(nr, 1))
284 CALL atom_density(lap(:, 1), atom%orbitals%pmat, atom%basis, atom%state%maxl_occ, &
285 typ="LAP", rr=atom%basis%grid%rad2)
286 IF (nlcc) THEN
287 CALL atom_core_density(lap(:, 1), atom%potential, typ="LAP", rr=atom%basis%grid%rad)
288 END IF
289 END IF
290
291 CALL fill_rho_set(rho_set, nspins, needs, rho, drho, tau, lap, nr)
292
293 CALL xc_dset_zero_all(deriv_set)
294
295 deriv_order = 1
296 CALL xc_functionals_eval(xc_fun_section, lsd=lsd, rho_set=rho_set, deriv_set=deriv_set, &
297 deriv_order=deriv_order)
298
299 ! Integration to get the matrix elements and energy
300 deriv => xc_dset_get_derivative(deriv_set, [INTEGER::], allocate_deriv=.false.)
301 CALL xc_derivative_get(deriv, deriv_data=xcpot)
302 atom%energy%exc = fourpi*integrate_grid(xcpot(:, 1, 1), atom%basis%grid)
303
304 ! dump grid density and xcpot (xc energy?)
305 IF (.false.) THEN
306 CALL open_file(unit_number=unit_nr, file_status="UNKNOWN", file_action="WRITE", file_name="atom.dat")
307 DO i = 1, SIZE(xcpot(:, 1, 1))
308 WRITE (unit_nr, *) atom%basis%grid%rad(i), rho(i, 1), xcpot(i, 1, 1)
309 END DO
310 CALL close_file(unit_nr)
311 END IF
312
313 IF (needs%rho) THEN
314 deriv => xc_dset_get_derivative(deriv_set, [deriv_rho], allocate_deriv=.false.)
315 CALL xc_derivative_get(deriv, deriv_data=xcpot)
316 CALL numpot_matrix(xcmat%op, xcpot(:, 1, 1), atom%basis, 0)
317 IF (.false.) THEN
318 CALL open_file(unit_number=unit_nr, file_status="UNKNOWN", file_action="WRITE", file_name="atompot.dat")
319 DO i = 1, SIZE(xcpot(:, 1, 1))
320 WRITE (unit_nr, *) atom%basis%grid%rad(i), rho(i, 1), xcpot(i, 1, 1)
321 END DO
322 CALL close_file(unit_nr)
323 END IF
324 DEALLOCATE (rho)
325 END IF
326 IF (needs%norm_drho) THEN
327 deriv => xc_dset_get_derivative(deriv_set, [deriv_norm_drho], allocate_deriv=.false.)
328 CALL xc_derivative_get(deriv, deriv_data=xcpot)
329 CALL numpot_matrix(xcmat%op, xcpot(:, 1, 1), atom%basis, 1)
330 IF (.false.) THEN
331 CALL open_file(unit_number=unit_nr, file_status="UNKNOWN", file_action="WRITE", file_name="atomdpot.dat")
332 DO i = 1, SIZE(xcpot(:, 1, 1))
333 WRITE (unit_nr, *) atom%basis%grid%rad(i), drho(i, 1), xcpot(i, 1, 1)
334 END DO
335 CALL close_file(unit_nr)
336 END IF
337 DEALLOCATE (drho)
338 END IF
339 IF (needs%tau) THEN
340 deriv => xc_dset_get_derivative(deriv_set, [deriv_tau], allocate_deriv=.false.)
341 CALL xc_derivative_get(deriv, deriv_data=xcpot)
342 n1 = SIZE(xcmat%op, 1)
343 n2 = SIZE(xcmat%op, 2)
344 n3 = SIZE(xcmat%op, 3)
345 ALLOCATE (taumat(n1, n2, 0:n3 - 1))
346 taumat = 0._dp
347
348 xcpot(:, 1, 1) = 0.5_dp*xcpot(:, 1, 1)
349 CALL numpot_matrix(xcmat%op, xcpot(:, 1, 1), atom%basis, 2)
350 xcpot(:, 1, 1) = xcpot(:, 1, 1)/atom%basis%grid%rad2(:)
351 CALL numpot_matrix(taumat, xcpot(:, 1, 1), atom%basis, 0)
352 DO l = 0, 3
353 xcmat%op(:, :, l) = xcmat%op(:, :, l) + real(l*(l + 1), dp)*taumat(:, :, l)
354 END DO
355
356 DEALLOCATE (tau)
357 DEALLOCATE (taumat)
358 END IF
359
360 ! assume lap is not really needed
361 IF (needs%laplace_rho) THEN
362 DEALLOCATE (lap)
363 END IF
364
365 ! Release the xc structure used to store the xc derivatives
366 CALL xc_dset_release(deriv_set)
367 CALL xc_rho_set_release(rho_set)
368
369 END IF !xc_none
370
371 ELSE
372
373 ! we don't have an xc_section, use a default setup
374 nr = atom%basis%grid%nr
375 ALLOCATE (rho(nr, 1), exc(nr), vxc(nr))
376
377 CALL atom_density(rho(:, 1), atom%orbitals%pmat, atom%basis, atom%state%maxl_occ, typ="RHO")
378 IF (nlcc) THEN
379 CALL atom_core_density(rho(:, 1), atom%potential, typ="RHO", rr=atom%basis%grid%rad)
380 END IF
381 CALL lda_pade(rho(:, 1), exc, vxc)
382
383 atom%energy%exc = fourpi*integrate_grid(exc, atom%basis%grid)
384 CALL numpot_matrix(xcmat%op, vxc, atom%basis, 0)
385
386 DEALLOCATE (rho, exc, vxc)
387
388 END IF
389
390 CALL timestop(handle)
391
392 END SUBROUTINE calculate_atom_vxc_lda
393
394! **************************************************************************************************
395!> \brief Calculate atomic exchange-correlation potential in spin-polarized case.
396!> \param xcmata spin-alpha exchange-correlation potential matrix
397!> \param xcmatb spin-beta exchange-correlation potential matrix
398!> \param atom information about the atomic kind
399!> \param xc_section XC input section
400!> \par History
401!> * 02.2010 created [Juerg Hutter]
402! **************************************************************************************************
403 SUBROUTINE calculate_atom_vxc_lsd(xcmata, xcmatb, atom, xc_section)
404 TYPE(opmat_type), POINTER :: xcmata, xcmatb
405 TYPE(atom_type), INTENT(INOUT) :: atom
406 TYPE(section_vals_type), POINTER :: xc_section
407
408 CHARACTER(LEN=*), PARAMETER :: routinen = 'calculate_atom_vxc_lsd'
409
410 INTEGER :: deriv_order, handle, l, myfun, n1, n2, &
411 n3, nr, nspins
412 INTEGER, DIMENSION(2, 3) :: bounds
413 LOGICAL :: lsd, nlcc
414 REAL(kind=dp) :: density_cut, gradient_cut, tau_cut
415 REAL(kind=dp), DIMENSION(:), POINTER :: exc, vxca, vxcb, xfun
416 REAL(kind=dp), DIMENSION(:, :), POINTER :: drho, lap, rho, tau
417 REAL(kind=dp), DIMENSION(:, :, :), POINTER :: taumat, xcpot
418 TYPE(section_vals_type), POINTER :: xc_fun_section
419 TYPE(xc_derivative_set_type) :: deriv_set
420 TYPE(xc_derivative_type), POINTER :: deriv
421 TYPE(xc_rho_cflags_type) :: needs
422 TYPE(xc_rho_set_type) :: rho_set
423
424 CALL timeset(routinen, handle)
425
426 nlcc = .false.
427 IF (atom%potential%ppot_type == gth_pseudo) THEN
428 nlcc = atom%potential%gth_pot%nlcc
429 ELSE IF (atom%potential%ppot_type == upf_pseudo) THEN
430 nlcc = atom%potential%upf_pot%core_correction
431 cpassert(.NOT. nlcc)
432 ELSE IF (atom%potential%ppot_type == sgp_pseudo) THEN
433 nlcc = atom%potential%sgp_pot%has_nlcc
434 cpassert(.NOT. nlcc)
435 END IF
436
437 IF (ASSOCIATED(xc_section)) THEN
438 IF (nlcc) THEN
439 nr = atom%basis%grid%nr
440 ALLOCATE (xfun(nr))
441 END IF
442
443 xc_fun_section => section_vals_get_subs_vals(xc_section, "XC_FUNCTIONAL")
444 CALL section_vals_val_get(xc_fun_section, "_SECTION_PARAMETERS_", i_val=myfun)
445
446 IF (myfun == xc_none) THEN
447 atom%energy%exc = 0._dp
448 ELSE
449 CALL section_vals_val_get(xc_section, "DENSITY_CUTOFF", r_val=density_cut)
450 CALL section_vals_val_get(xc_section, "GRADIENT_CUTOFF", r_val=gradient_cut)
451 CALL section_vals_val_get(xc_section, "TAU_CUTOFF", r_val=tau_cut)
452
453 lsd = .true.
454 nspins = 2
455 needs = xc_functionals_get_needs(xc_fun_section, lsd=lsd, calc_potential=.false.)
456
457 ! Prepare the structures needed to calculate and store the xc derivatives
458
459 ! Array dimension: here anly one dimensional arrays are used,
460 ! i.e. only the first column of deriv_data is read.
461 ! The other to dimensions are set to size equal 1
462 nr = atom%basis%grid%nr
463 bounds(1:2, 1:3) = 1
464 bounds(2, 1) = nr
465
466 ! create a place where to put the derivatives
467 CALL xc_dset_create(deriv_set, local_bounds=bounds)
468 ! create the place where to store the argument for the functionals
469 CALL xc_rho_set_create(rho_set, bounds, rho_cutoff=density_cut, &
470 drho_cutoff=gradient_cut, tau_cutoff=tau_cut)
471 ! allocate the required 3d arrays where to store rho and drho
472 CALL xc_rho_set_atom_update(rho_set, needs, nspins, bounds)
473
474 NULLIFY (rho, drho, tau)
475 IF (needs%rho_spin) THEN
476 ALLOCATE (rho(nr, 2))
477 CALL atom_density(rho(:, 1), atom%orbitals%pmata, atom%basis, atom%state%maxl_occ, typ="RHO")
478 CALL atom_density(rho(:, 2), atom%orbitals%pmatb, atom%basis, atom%state%maxl_occ, typ="RHO")
479 IF (nlcc) THEN
480 xfun = 0.0_dp
481 CALL atom_core_density(xfun(:), atom%potential, typ="RHO", rr=atom%basis%grid%rad)
482 rho(:, 1) = rho(:, 1) + 0.5_dp*xfun(:)
483 rho(:, 2) = rho(:, 2) + 0.5_dp*xfun(:)
484 END IF
485 END IF
486 IF (needs%norm_drho_spin) THEN
487 ALLOCATE (drho(nr, 2))
488 CALL atom_density(drho(:, 1), atom%orbitals%pmata, atom%basis, atom%state%maxl_occ, typ="DER")
489 CALL atom_density(drho(:, 2), atom%orbitals%pmatb, atom%basis, atom%state%maxl_occ, typ="DER")
490 IF (nlcc) THEN
491 xfun = 0.0_dp
492 CALL atom_core_density(xfun(:), atom%potential, typ="DER", rr=atom%basis%grid%rad)
493 drho(:, 1) = drho(:, 1) + 0.5_dp*xfun(:)
494 drho(:, 2) = drho(:, 2) + 0.5_dp*xfun(:)
495 END IF
496 END IF
497 IF (needs%tau_spin) THEN
498 ALLOCATE (tau(nr, 2))
499 CALL atom_density(tau(:, 1), atom%orbitals%pmata, atom%basis, atom%state%maxl_occ, &
500 typ="KIN", rr=atom%basis%grid%rad2)
501 CALL atom_density(tau(:, 2), atom%orbitals%pmatb, atom%basis, atom%state%maxl_occ, &
502 typ="KIN", rr=atom%basis%grid%rad2)
503 END IF
504 IF (needs%laplace_rho_spin) THEN
505 ALLOCATE (lap(nr, 2))
506 CALL atom_density(lap(:, 1), atom%orbitals%pmata, atom%basis, atom%state%maxl_occ, &
507 typ="LAP", rr=atom%basis%grid%rad2)
508 CALL atom_density(lap(:, 2), atom%orbitals%pmatb, atom%basis, atom%state%maxl_occ, &
509 typ="LAP", rr=atom%basis%grid%rad2)
510 IF (nlcc) THEN
511 xfun = 0.0_dp
512 CALL atom_core_density(xfun(:), atom%potential, typ="LAP", rr=atom%basis%grid%rad)
513 lap(:, 1) = lap(:, 1) + 0.5_dp*xfun(:)
514 lap(:, 2) = lap(:, 2) + 0.5_dp*xfun(:)
515 END IF
516 END IF
517
518 CALL fill_rho_set(rho_set, nspins, needs, rho, drho, tau, lap, nr)
519
520 CALL xc_dset_zero_all(deriv_set)
521
522 deriv_order = 1
523 CALL xc_functionals_eval(xc_fun_section, lsd=lsd, rho_set=rho_set, deriv_set=deriv_set, &
524 deriv_order=deriv_order)
525
526 ! Integration to get the matrix elements and energy
527 deriv => xc_dset_get_derivative(deriv_set, [INTEGER::], allocate_deriv=.false.)
528 CALL xc_derivative_get(deriv, deriv_data=xcpot)
529 atom%energy%exc = fourpi*integrate_grid(xcpot(:, 1, 1), atom%basis%grid)
530
531 IF (needs%rho_spin) THEN
532 deriv => xc_dset_get_derivative(deriv_set, [deriv_rhoa], allocate_deriv=.false.)
533 CALL xc_derivative_get(deriv, deriv_data=xcpot)
534 CALL numpot_matrix(xcmata%op, xcpot(:, 1, 1), atom%basis, 0)
535 deriv => xc_dset_get_derivative(deriv_set, [deriv_rhob], allocate_deriv=.false.)
536 CALL xc_derivative_get(deriv, deriv_data=xcpot)
537 CALL numpot_matrix(xcmatb%op, xcpot(:, 1, 1), atom%basis, 0)
538 DEALLOCATE (rho)
539 END IF
540 IF (needs%norm_drho_spin) THEN
541 ! drhoa
542 NULLIFY (deriv)
543 deriv => xc_dset_get_derivative(deriv_set, [deriv_norm_drhoa], allocate_deriv=.false.)
544 CALL xc_derivative_get(deriv, deriv_data=xcpot)
545 CALL numpot_matrix(xcmata%op, xcpot(:, 1, 1), atom%basis, 1)
546 ! drhob
547 NULLIFY (deriv)
548 deriv => xc_dset_get_derivative(deriv_set, [deriv_norm_drhob], allocate_deriv=.false.)
549 CALL xc_derivative_get(deriv, deriv_data=xcpot)
550 CALL numpot_matrix(xcmatb%op, xcpot(:, 1, 1), atom%basis, 1)
551 ! Cross Terms
552 NULLIFY (deriv)
553 deriv => xc_dset_get_derivative(deriv_set, [deriv_norm_drho])
554 IF (ASSOCIATED(deriv)) THEN
555 CALL xc_derivative_get(deriv, deriv_data=xcpot)
556 CALL numpot_matrix(xcmata%op, xcpot(:, 1, 1), atom%basis, 1)
557 CALL numpot_matrix(xcmatb%op, xcpot(:, 1, 1), atom%basis, 1)
558 END IF
559 DEALLOCATE (drho)
560 END IF
561 IF (needs%tau_spin) THEN
562 n1 = SIZE(xcmata%op, 1)
563 n2 = SIZE(xcmata%op, 2)
564 n3 = SIZE(xcmata%op, 3)
565 ALLOCATE (taumat(n1, n2, 0:n3 - 1))
566
567 deriv => xc_dset_get_derivative(deriv_set, [deriv_tau_a], allocate_deriv=.false.)
568 CALL xc_derivative_get(deriv, deriv_data=xcpot)
569 taumat = 0._dp
570 xcpot(:, 1, 1) = 0.5_dp*xcpot(:, 1, 1)
571 CALL numpot_matrix(xcmata%op, xcpot(:, 1, 1), atom%basis, 2)
572 xcpot(:, 1, 1) = xcpot(:, 1, 1)/atom%basis%grid%rad2(:)
573 CALL numpot_matrix(taumat, xcpot(:, 1, 1), atom%basis, 0)
574 DO l = 0, 3
575 xcmata%op(:, :, l) = xcmata%op(:, :, l) + real(l*(l + 1), dp)*taumat(:, :, l)
576 END DO
577
578 deriv => xc_dset_get_derivative(deriv_set, [deriv_tau_b], allocate_deriv=.false.)
579 CALL xc_derivative_get(deriv, deriv_data=xcpot)
580 taumat = 0._dp
581 xcpot(:, 1, 1) = 0.5_dp*xcpot(:, 1, 1)
582 CALL numpot_matrix(xcmatb%op, xcpot(:, 1, 1), atom%basis, 2)
583 xcpot(:, 1, 1) = xcpot(:, 1, 1)/atom%basis%grid%rad2(:)
584 CALL numpot_matrix(taumat, xcpot(:, 1, 1), atom%basis, 0)
585 DO l = 0, 3
586 xcmatb%op(:, :, l) = xcmatb%op(:, :, l) + real(l*(l + 1), dp)*taumat(:, :, l)
587 END DO
588
589 DEALLOCATE (tau)
590 DEALLOCATE (taumat)
591 END IF
592
593 ! assume lap is not really needed
594 IF (needs%laplace_rho_spin) THEN
595 DEALLOCATE (lap)
596 END IF
597
598 ! Release the xc structure used to store the xc derivatives
599 CALL xc_dset_release(deriv_set)
600 CALL xc_rho_set_release(rho_set)
601
602 END IF !xc_none
603
604 IF (nlcc) DEALLOCATE (xfun)
605
606 ELSE
607
608 ! we don't have an xc_section, use a default setup
609 nr = atom%basis%grid%nr
610 ALLOCATE (rho(nr, 2), exc(nr), vxca(nr), vxcb(nr))
611 IF (nlcc) ALLOCATE (xfun(nr))
612
613 CALL atom_density(rho(:, 1), atom%orbitals%pmata, atom%basis, atom%state%maxl_occ, typ="RHO")
614 CALL atom_density(rho(:, 2), atom%orbitals%pmatb, atom%basis, atom%state%maxl_occ, typ="RHO")
615 IF (nlcc) THEN
616 xfun(:) = 0.0_dp
617 CALL atom_core_density(xfun(:), atom%potential, typ="RHO", rr=atom%basis%grid%rad)
618 rho(:, 1) = rho(:, 1) + 0.5_dp*xfun(:)
619 rho(:, 2) = rho(:, 2) + 0.5_dp*xfun(:)
620 END IF
621 CALL lsd_pade(rho(:, 1), rho(:, 2), exc, vxca, vxcb)
622
623 atom%energy%exc = fourpi*integrate_grid(exc, atom%basis%grid)
624 CALL numpot_matrix(xcmata%op, vxca, atom%basis, 0)
625 CALL numpot_matrix(xcmatb%op, vxcb, atom%basis, 0)
626
627 DEALLOCATE (rho, exc, vxca, vxcb)
628 IF (nlcc) DEALLOCATE (xfun)
629
630 END IF
631
632 CALL timestop(handle)
633
634 END SUBROUTINE calculate_atom_vxc_lsd
635
636! **************************************************************************************************
637!> \brief Alternative subroutine that calculates atomic exchange-correlation potential
638!> in spin-restricted case.
639!> \param basis atomic basis set
640!> \param pmat density matrix
641!> \param maxl maximum angular momentum
642!> \param xc_section XC input section
643!> \param fexc exchange-correlation energy
644!> \param xcmat exchange-correlation potential matrix
645!> \par History
646!> * 07.2016 ADMM analysis [Juerg Hutter]
647! **************************************************************************************************
648 SUBROUTINE atom_vxc_lda(basis, pmat, maxl, xc_section, fexc, xcmat)
649 TYPE(atom_basis_type), POINTER :: basis
650 REAL(kind=dp), DIMENSION(:, :, :), INTENT(IN) :: pmat
651 INTEGER, INTENT(IN) :: maxl
652 TYPE(section_vals_type), POINTER :: xc_section
653 REAL(kind=dp), INTENT(OUT) :: fexc
654 REAL(kind=dp), DIMENSION(:, :, :), INTENT(INOUT) :: xcmat
655
656 CHARACTER(LEN=*), PARAMETER :: routinen = 'atom_vxc_lda'
657
658 INTEGER :: deriv_order, handle, l, myfun, n1, n2, &
659 n3, nr, nspins
660 INTEGER, DIMENSION(2, 3) :: bounds
661 LOGICAL :: lsd
662 REAL(kind=dp) :: density_cut, gradient_cut, tau_cut
663 REAL(kind=dp), DIMENSION(:, :), POINTER :: drho, lap, rho, tau
664 REAL(kind=dp), DIMENSION(:, :, :), POINTER :: taumat, xcpot
665 TYPE(section_vals_type), POINTER :: xc_fun_section
666 TYPE(xc_derivative_set_type) :: deriv_set
667 TYPE(xc_derivative_type), POINTER :: deriv
668 TYPE(xc_rho_cflags_type) :: needs
669 TYPE(xc_rho_set_type) :: rho_set
670
671 CALL timeset(routinen, handle)
672
673 cpassert(ASSOCIATED(xc_section))
674
675 xc_fun_section => section_vals_get_subs_vals(xc_section, "XC_FUNCTIONAL")
676 CALL section_vals_val_get(xc_fun_section, "_SECTION_PARAMETERS_", i_val=myfun)
677
678 IF (myfun == xc_none) THEN
679 fexc = 0._dp
680 ELSE
681 CALL section_vals_val_get(xc_section, "DENSITY_CUTOFF", r_val=density_cut)
682 CALL section_vals_val_get(xc_section, "GRADIENT_CUTOFF", r_val=gradient_cut)
683 CALL section_vals_val_get(xc_section, "TAU_CUTOFF", r_val=tau_cut)
684
685 lsd = .false.
686 nspins = 1
687 needs = xc_functionals_get_needs(xc_fun_section, lsd=lsd, calc_potential=.false.)
688
689 ! Prepare the structures needed to calculate and store the xc derivatives
690
691 ! Array dimension: here anly one dimensional arrays are used,
692 ! i.e. only the first column of deriv_data is read.
693 ! The other to dimensions are set to size equal 1
694 nr = basis%grid%nr
695 bounds(1:2, 1:3) = 1
696 bounds(2, 1) = nr
697
698 ! create a place where to put the derivatives
699 CALL xc_dset_create(deriv_set, local_bounds=bounds)
700 ! create the place where to store the argument for the functionals
701 CALL xc_rho_set_create(rho_set, bounds, rho_cutoff=density_cut, &
702 drho_cutoff=gradient_cut, tau_cutoff=tau_cut)
703 ! allocate the required 3d arrays where to store rho and drho
704 CALL xc_rho_set_atom_update(rho_set, needs, nspins, bounds)
705
706 NULLIFY (rho, drho, tau)
707 IF (needs%rho) THEN
708 ALLOCATE (rho(nr, 1))
709 CALL atom_density(rho(:, 1), pmat, basis, maxl, typ="RHO")
710 END IF
711 IF (needs%norm_drho) THEN
712 ALLOCATE (drho(nr, 1))
713 CALL atom_density(drho(:, 1), pmat, basis, maxl, typ="DER")
714 END IF
715 IF (needs%tau) THEN
716 ALLOCATE (tau(nr, 1))
717 CALL atom_density(tau(:, 1), pmat, basis, maxl, typ="KIN", rr=basis%grid%rad2)
718 END IF
719 IF (needs%laplace_rho) THEN
720 ALLOCATE (lap(nr, 1))
721 CALL atom_density(lap(:, 1), pmat, basis, maxl, typ="LAP", rr=basis%grid%rad2)
722 END IF
723
724 CALL fill_rho_set(rho_set, nspins, needs, rho, drho, tau, lap, nr)
725
726 CALL xc_dset_zero_all(deriv_set)
727
728 deriv_order = 1
729 CALL xc_functionals_eval(xc_fun_section, lsd=lsd, rho_set=rho_set, deriv_set=deriv_set, &
730 deriv_order=deriv_order)
731
732 ! Integration to get the matrix elements and energy
733 deriv => xc_dset_get_derivative(deriv_set, [INTEGER::], allocate_deriv=.false.)
734 CALL xc_derivative_get(deriv, deriv_data=xcpot)
735 fexc = fourpi*integrate_grid(xcpot(:, 1, 1), basis%grid)
736
737 IF (needs%rho) THEN
738 deriv => xc_dset_get_derivative(deriv_set, [deriv_rho], allocate_deriv=.false.)
739 CALL xc_derivative_get(deriv, deriv_data=xcpot)
740 CALL numpot_matrix(xcmat, xcpot(:, 1, 1), basis, 0)
741 DEALLOCATE (rho)
742 END IF
743 IF (needs%norm_drho) THEN
744 deriv => xc_dset_get_derivative(deriv_set, [deriv_norm_drho], allocate_deriv=.false.)
745 CALL xc_derivative_get(deriv, deriv_data=xcpot)
746 CALL numpot_matrix(xcmat, xcpot(:, 1, 1), basis, 1)
747 DEALLOCATE (drho)
748 END IF
749 IF (needs%tau) THEN
750 deriv => xc_dset_get_derivative(deriv_set, [deriv_tau], allocate_deriv=.false.)
751 CALL xc_derivative_get(deriv, deriv_data=xcpot)
752 n1 = SIZE(xcmat, 1)
753 n2 = SIZE(xcmat, 2)
754 n3 = SIZE(xcmat, 3)
755 ALLOCATE (taumat(n1, n2, 0:n3 - 1))
756 taumat = 0._dp
757
758 xcpot(:, 1, 1) = 0.5_dp*xcpot(:, 1, 1)
759 CALL numpot_matrix(xcmat, xcpot(:, 1, 1), basis, 2)
760 xcpot(:, 1, 1) = xcpot(:, 1, 1)/basis%grid%rad2(:)
761 CALL numpot_matrix(taumat, xcpot(:, 1, 1), basis, 0)
762 DO l = 0, 3
763 xcmat(:, :, l) = xcmat(:, :, l) + real(l*(l + 1), dp)*taumat(:, :, l)
764 END DO
765
766 DEALLOCATE (tau)
767 DEALLOCATE (taumat)
768 END IF
769
770 ! we assume lap is not really needed
771 IF (needs%laplace_rho) THEN
772 DEALLOCATE (lap)
773 END IF
774
775 ! Release the xc structure used to store the xc derivatives
776 CALL xc_dset_release(deriv_set)
777 CALL xc_rho_set_release(rho_set)
778
779 END IF !xc_none
780
781 CALL timestop(handle)
782
783 END SUBROUTINE atom_vxc_lda
784
785! **************************************************************************************************
786!> \brief Alternative subroutine that calculates atomic exchange-correlation potential
787!> in spin-polarized case.
788!> \param basis atomic basis set
789!> \param pmata spin-alpha density matrix
790!> \param pmatb spin-beta density matrix
791!> \param maxl maximum angular momentum
792!> \param xc_section XC input section
793!> \param fexc exchange-correlation energy
794!> \param xcmata spin-alpha exchange-correlation potential matrix
795!> \param xcmatb spin-beta exchange-correlation potential matrix
796!> \par History
797!> * 07.2016 ADMM analysis [Juerg Hutter]
798! **************************************************************************************************
799 SUBROUTINE atom_vxc_lsd(basis, pmata, pmatb, maxl, xc_section, fexc, xcmata, xcmatb)
800 TYPE(atom_basis_type), POINTER :: basis
801 REAL(kind=dp), DIMENSION(:, :, :), INTENT(IN) :: pmata, pmatb
802 INTEGER, INTENT(IN) :: maxl
803 TYPE(section_vals_type), POINTER :: xc_section
804 REAL(kind=dp), INTENT(OUT) :: fexc
805 REAL(kind=dp), DIMENSION(:, :, :), INTENT(INOUT) :: xcmata, xcmatb
806
807 CHARACTER(LEN=*), PARAMETER :: routinen = 'atom_vxc_lsd'
808
809 INTEGER :: deriv_order, handle, l, myfun, n1, n2, &
810 n3, nr, nspins
811 INTEGER, DIMENSION(2, 3) :: bounds
812 LOGICAL :: lsd
813 REAL(kind=dp) :: density_cut, gradient_cut, tau_cut
814 REAL(kind=dp), DIMENSION(:, :), POINTER :: drho, lap, rho, tau
815 REAL(kind=dp), DIMENSION(:, :, :), POINTER :: taumat, xcpot
816 TYPE(section_vals_type), POINTER :: xc_fun_section
817 TYPE(xc_derivative_set_type) :: deriv_set
818 TYPE(xc_derivative_type), POINTER :: deriv
819 TYPE(xc_rho_cflags_type) :: needs
820 TYPE(xc_rho_set_type) :: rho_set
821
822 CALL timeset(routinen, handle)
823
824 cpassert(ASSOCIATED(xc_section))
825
826 xc_fun_section => section_vals_get_subs_vals(xc_section, "XC_FUNCTIONAL")
827 CALL section_vals_val_get(xc_fun_section, "_SECTION_PARAMETERS_", i_val=myfun)
828
829 IF (myfun == xc_none) THEN
830 fexc = 0._dp
831 ELSE
832 CALL section_vals_val_get(xc_section, "DENSITY_CUTOFF", r_val=density_cut)
833 CALL section_vals_val_get(xc_section, "GRADIENT_CUTOFF", r_val=gradient_cut)
834 CALL section_vals_val_get(xc_section, "TAU_CUTOFF", r_val=tau_cut)
835
836 lsd = .true.
837 nspins = 2
838 needs = xc_functionals_get_needs(xc_fun_section, lsd=lsd, calc_potential=.false.)
839
840 ! Prepare the structures needed to calculate and store the xc derivatives
841
842 ! Array dimension: here anly one dimensional arrays are used,
843 ! i.e. only the first column of deriv_data is read.
844 ! The other to dimensions are set to size equal 1
845 nr = basis%grid%nr
846 bounds(1:2, 1:3) = 1
847 bounds(2, 1) = nr
848
849 ! create a place where to put the derivatives
850 CALL xc_dset_create(deriv_set, local_bounds=bounds)
851 ! create the place where to store the argument for the functionals
852 CALL xc_rho_set_create(rho_set, bounds, rho_cutoff=density_cut, &
853 drho_cutoff=gradient_cut, tau_cutoff=tau_cut)
854 ! allocate the required 3d arrays where to store rho and drho
855 CALL xc_rho_set_atom_update(rho_set, needs, nspins, bounds)
856
857 NULLIFY (rho, drho, tau)
858 IF (needs%rho_spin) THEN
859 ALLOCATE (rho(nr, 2))
860 CALL atom_density(rho(:, 1), pmata, basis, maxl, typ="RHO")
861 CALL atom_density(rho(:, 2), pmatb, basis, maxl, typ="RHO")
862 END IF
863 IF (needs%norm_drho_spin) THEN
864 ALLOCATE (drho(nr, 2))
865 CALL atom_density(drho(:, 1), pmata, basis, maxl, typ="DER")
866 CALL atom_density(drho(:, 2), pmatb, basis, maxl, typ="DER")
867 END IF
868 IF (needs%tau_spin) THEN
869 ALLOCATE (tau(nr, 2))
870 CALL atom_density(tau(:, 1), pmata, basis, maxl, typ="KIN", rr=basis%grid%rad2)
871 CALL atom_density(tau(:, 2), pmatb, basis, maxl, typ="KIN", rr=basis%grid%rad2)
872 END IF
873 IF (needs%laplace_rho_spin) THEN
874 ALLOCATE (lap(nr, 2))
875 CALL atom_density(lap(:, 1), pmata, basis, maxl, typ="LAP", rr=basis%grid%rad2)
876 CALL atom_density(lap(:, 2), pmatb, basis, maxl, typ="LAP", rr=basis%grid%rad2)
877 END IF
878
879 CALL fill_rho_set(rho_set, nspins, needs, rho, drho, tau, lap, nr)
880
881 CALL xc_dset_zero_all(deriv_set)
882
883 deriv_order = 1
884 CALL xc_functionals_eval(xc_fun_section, lsd=lsd, rho_set=rho_set, deriv_set=deriv_set, &
885 deriv_order=deriv_order)
886
887 ! Integration to get the matrix elements and energy
888 deriv => xc_dset_get_derivative(deriv_set, [INTEGER::], allocate_deriv=.false.)
889 CALL xc_derivative_get(deriv, deriv_data=xcpot)
890 fexc = fourpi*integrate_grid(xcpot(:, 1, 1), basis%grid)
891
892 IF (needs%rho_spin) THEN
893 deriv => xc_dset_get_derivative(deriv_set, [deriv_rhoa], allocate_deriv=.false.)
894 CALL xc_derivative_get(deriv, deriv_data=xcpot)
895 CALL numpot_matrix(xcmata, xcpot(:, 1, 1), basis, 0)
896 deriv => xc_dset_get_derivative(deriv_set, [deriv_rhob], allocate_deriv=.false.)
897 CALL xc_derivative_get(deriv, deriv_data=xcpot)
898 CALL numpot_matrix(xcmatb, xcpot(:, 1, 1), basis, 0)
899 DEALLOCATE (rho)
900 END IF
901 IF (needs%norm_drho_spin) THEN
902 ! drhoa
903 NULLIFY (deriv)
904 deriv => xc_dset_get_derivative(deriv_set, [deriv_norm_drhoa], allocate_deriv=.false.)
905 CALL xc_derivative_get(deriv, deriv_data=xcpot)
906 CALL numpot_matrix(xcmata, xcpot(:, 1, 1), basis, 1)
907 ! drhob
908 NULLIFY (deriv)
909 deriv => xc_dset_get_derivative(deriv_set, [deriv_norm_drhob], allocate_deriv=.false.)
910 CALL xc_derivative_get(deriv, deriv_data=xcpot)
911 CALL numpot_matrix(xcmatb, xcpot(:, 1, 1), basis, 1)
912 ! Cross Terms
913 NULLIFY (deriv)
914 deriv => xc_dset_get_derivative(deriv_set, [deriv_norm_drho])
915 IF (ASSOCIATED(deriv)) THEN
916 CALL xc_derivative_get(deriv, deriv_data=xcpot)
917 CALL numpot_matrix(xcmata, xcpot(:, 1, 1), basis, 1)
918 CALL numpot_matrix(xcmatb, xcpot(:, 1, 1), basis, 1)
919 END IF
920 DEALLOCATE (drho)
921 END IF
922 IF (needs%tau_spin) THEN
923 n1 = SIZE(xcmata, 1)
924 n2 = SIZE(xcmata, 2)
925 n3 = SIZE(xcmata, 3)
926 ALLOCATE (taumat(n1, n2, 0:n3 - 1))
927
928 deriv => xc_dset_get_derivative(deriv_set, [deriv_tau_a], allocate_deriv=.false.)
929 CALL xc_derivative_get(deriv, deriv_data=xcpot)
930 taumat = 0._dp
931 xcpot(:, 1, 1) = 0.5_dp*xcpot(:, 1, 1)
932 CALL numpot_matrix(xcmata, xcpot(:, 1, 1), basis, 2)
933 xcpot(:, 1, 1) = xcpot(:, 1, 1)/basis%grid%rad2(:)
934 CALL numpot_matrix(taumat, xcpot(:, 1, 1), basis, 0)
935 DO l = 0, 3
936 xcmata(:, :, l) = xcmata(:, :, l) + real(l*(l + 1), dp)*taumat(:, :, l)
937 END DO
938
939 deriv => xc_dset_get_derivative(deriv_set, [deriv_tau_b], allocate_deriv=.false.)
940 CALL xc_derivative_get(deriv, deriv_data=xcpot)
941 taumat = 0._dp
942 xcpot(:, 1, 1) = 0.5_dp*xcpot(:, 1, 1)
943 CALL numpot_matrix(xcmatb, xcpot(:, 1, 1), basis, 2)
944 xcpot(:, 1, 1) = xcpot(:, 1, 1)/basis%grid%rad2(:)
945 CALL numpot_matrix(taumat, xcpot(:, 1, 1), basis, 0)
946 DO l = 0, 3
947 xcmatb(:, :, l) = xcmatb(:, :, l) + real(l*(l + 1), dp)*taumat(:, :, l)
948 END DO
949
950 DEALLOCATE (tau)
951 DEALLOCATE (taumat)
952 END IF
953
954 ! Release the xc structure used to store the xc derivatives
955 CALL xc_dset_release(deriv_set)
956 CALL xc_rho_set_release(rho_set)
957
958 END IF !xc_none
959
960 CALL timestop(handle)
961
962 END SUBROUTINE atom_vxc_lsd
963
964! **************************************************************************************************
965!> \brief Estimate the ADMM exchange energy correction term using a model exchange functional.
966!> \param basis0 reference atomic basis set
967!> \param pmat0 reference density matrix
968!> \param basis1 ADMM basis set
969!> \param pmat1 ADMM density matrix
970!> \param maxl maxumum angular momentum
971!> \param functional name of the model exchange functional:
972!> "LINX" -- linear extrapolation of the Slater exchange potential [?]
973!> \param dfexc exchange-correlation energy difference
974!> \param linxpar LINX parameters
975!> \par History
976!> * 07.2016 ADMM analysis [Juerg Hutter]
977! **************************************************************************************************
978 SUBROUTINE atom_dpot_lda(basis0, pmat0, basis1, pmat1, maxl, functional, dfexc, linxpar)
979 TYPE(atom_basis_type), POINTER :: basis0
980 REAL(kind=dp), DIMENSION(:, :, :), INTENT(IN) :: pmat0
981 TYPE(atom_basis_type), POINTER :: basis1
982 REAL(kind=dp), DIMENSION(:, :, :), INTENT(IN) :: pmat1
983 INTEGER, INTENT(IN) :: maxl
984 CHARACTER(LEN=*), INTENT(IN) :: functional
985 REAL(kind=dp), INTENT(OUT) :: dfexc
986 REAL(kind=dp), DIMENSION(:), INTENT(IN), OPTIONAL :: linxpar
987
988 CHARACTER(LEN=*), PARAMETER :: routinen = 'atom_dpot_lda'
989 REAL(kind=dp), PARAMETER :: alx = 0.20_dp, blx = -0.06_dp, &
990 fs = -0.75_dp*(3._dp/pi)**(1._dp/3._dp)
991
992 INTEGER :: handle, ir, nr
993 REAL(kind=dp) :: a, b, fx
994 REAL(kind=dp), ALLOCATABLE, DIMENSION(:) :: delta, drho0, drho1, pot0, pot1, rho0, &
995 rho1
996
997 CALL timeset(routinen, handle)
998
999 nr = basis0%grid%nr
1000
1001 ALLOCATE (rho0(nr), drho0(nr))
1002 CALL atom_density(rho0, pmat0, basis0, maxl, typ="RHO")
1003 CALL atom_density(drho0, pmat0, basis0, maxl, typ="DER")
1004 !
1005 ALLOCATE (rho1(nr), drho1(nr))
1006 CALL atom_density(rho1, pmat1, basis1, maxl, typ="RHO")
1007 ! CONSIDER TO REMOVE [?]: drho1 is calculated but it is not used.
1008 CALL atom_density(drho1, pmat1, basis1, maxl, typ="DER")
1009 !
1010 ALLOCATE (delta(nr))
1011 fx = 4.0_dp/3.0_dp
1012 delta(1:nr) = fs*(rho0(1:nr)**fx - rho1(1:nr)**fx)
1013
1014 SELECT CASE (trim(functional))
1015 CASE ("LINX")
1016 IF (PRESENT(linxpar)) THEN
1017 a = linxpar(1)
1018 b = linxpar(2)
1019 ELSE
1020 a = alx
1021 b = blx
1022 END IF
1023 ALLOCATE (pot0(nr), pot1(nr))
1024 DO ir = 1, nr
1025 IF (rho0(ir) > 1.e-12_dp) THEN
1026 pot0(ir) = 0.5_dp*drho0(ir)/(3._dp*pi*pi*rho0(ir)**fx)
1027 ELSE
1028 pot0(ir) = 0._dp
1029 END IF
1030 END DO
1031 pot1(1:nr) = 1._dp + (a*pot0(1:nr)**2)/(1._dp + b*pot0(1:nr)**2)
1032 pot1(1:nr) = pot1(1:nr)*delta(1:nr)
1033 dfexc = fourpi*integrate_grid(pot1(1:nr), basis0%grid)
1034 DEALLOCATE (pot0, pot1)
1035 CASE DEFAULT
1036 cpabort("Unknown functional in atom_dpot_lda")
1037 END SELECT
1038
1039 DEALLOCATE (rho0, rho1, drho0, drho1, delta)
1040
1041 CALL timestop(handle)
1042
1043 END SUBROUTINE atom_dpot_lda
1044
1045! **************************************************************************************************
1046!> \brief Initialise object that holds components needed to compute the exchange-correlation
1047!> functional on the atomic radial grid.
1048!> \param rho_set electron density object to initialise
1049!> \param nspins number of spin components
1050!> \param needs components needed to compute the exchange-correlation functional
1051!> \param rho electron density on the grid
1052!> \param drho first derivative of the electron density on the grid
1053!> \param tau kinetic energy density on the grid
1054!> \param lap second derivative of the electron density on the grid
1055!> \param na number of points on the atomic radial grid
1056!> \par History
1057!> * 02.2010 unrestricted KS and HF methods [Juerg Hutter]
1058!> * 08.2008 created as calculate_atom_vxc_lda() [Juerg Hutter]
1059! **************************************************************************************************
1060 SUBROUTINE fill_rho_set(rho_set, nspins, needs, rho, drho, tau, lap, na)
1061
1062 TYPE(xc_rho_set_type), INTENT(INOUT) :: rho_set
1063 INTEGER, INTENT(IN) :: nspins
1064 TYPE(xc_rho_cflags_type), INTENT(in) :: needs
1065 REAL(dp), DIMENSION(:, :), POINTER :: rho, drho, tau, lap
1066 INTEGER, INTENT(IN) :: na
1067
1068 REAL(kind=dp), PARAMETER :: f13 = (1.0_dp/3.0_dp)
1069
1070 INTEGER :: ia
1071
1072 SELECT CASE (nspins)
1073 CASE (1)
1074 cpassert(.NOT. needs%rho_spin)
1075 cpassert(.NOT. needs%drho_spin)
1076 cpassert(.NOT. needs%norm_drho_spin)
1077 cpassert(.NOT. needs%rho_spin_1_3)
1078 cpassert(.NOT. needs%tau_spin)
1079 cpassert(.NOT. needs%drho)
1080 ! Give rho to 1/3
1081 IF (needs%rho_1_3) THEN
1082 DO ia = 1, na
1083 rho_set%rho_1_3(ia, 1, 1) = max(rho(ia, 1), 0.0_dp)**f13
1084 END DO
1085 rho_set%owns%rho_1_3 = .true.
1086 rho_set%has%rho_1_3 = .true.
1087 END IF
1088 ! Give the density
1089 IF (needs%rho) THEN
1090 DO ia = 1, na
1091 rho_set%rho(ia, 1, 1) = rho(ia, 1)
1092 END DO
1093 rho_set%owns%rho = .true.
1094 rho_set%has%rho = .true.
1095 END IF
1096 ! Give the norm of the gradient of the density
1097 IF (needs%norm_drho) THEN
1098 DO ia = 1, na
1099 rho_set%norm_drho(ia, 1, 1) = drho(ia, 1)
1100 END DO
1101 rho_set%owns%norm_drho = .true.
1102 rho_set%has%norm_drho = .true.
1103 END IF
1104 CASE (2)
1105 cpassert(.NOT. needs%drho)
1106 cpassert(.NOT. needs%drho_spin)
1107 ! Give the total density
1108 IF (needs%rho) THEN
1109 DO ia = 1, na
1110 rho_set%rho(ia, 1, 1) = rho(ia, 1) + rho(ia, 2)
1111 END DO
1112 rho_set%owns%rho = .true.
1113 rho_set%has%rho = .true.
1114 END IF
1115 ! Give the norm of the total gradient of the density
1116 IF (needs%norm_drho) THEN
1117 DO ia = 1, na
1118 rho_set%norm_drho(ia, 1, 1) = drho(ia, 1) + drho(ia, 2)
1119 END DO
1120 rho_set%owns%norm_drho = .true.
1121 rho_set%has%norm_drho = .true.
1122 END IF
1123 ! Give rho_spin
1124 IF (needs%rho_spin) THEN
1125 DO ia = 1, na
1126 rho_set%rhoa(ia, 1, 1) = rho(ia, 1)
1127 rho_set%rhob(ia, 1, 1) = rho(ia, 2)
1128 END DO
1129 rho_set%owns%rho_spin = .true.
1130 rho_set%has%rho_spin = .true.
1131 END IF
1132 ! Give rho_spin to 1/3
1133 IF (needs%rho_spin_1_3) THEN
1134 DO ia = 1, na
1135 rho_set%rhoa_1_3(ia, 1, 1) = max(rho(ia, 1), 0.0_dp)**f13
1136 rho_set%rhob_1_3(ia, 1, 1) = max(rho(ia, 2), 0.0_dp)**f13
1137 END DO
1138 rho_set%owns%rho_1_3 = .true.
1139 rho_set%has%rho_1_3 = .true.
1140 END IF
1141 ! Give the norm of the gradient of rhoa and of rhob separatedly
1142 IF (needs%norm_drho_spin) THEN
1143 DO ia = 1, na
1144 rho_set%norm_drhoa(ia, 1, 1) = drho(ia, 1)
1145 rho_set%norm_drhob(ia, 1, 1) = drho(ia, 2)
1146 END DO
1147 rho_set%owns%norm_drho_spin = .true.
1148 rho_set%has%norm_drho_spin = .true.
1149 END IF
1150 !
1151 END SELECT
1152
1153 ! tau part
1154 IF (needs%tau) THEN
1155 IF (nspins == 2) THEN
1156 DO ia = 1, na
1157 rho_set%tau(ia, 1, 1) = tau(ia, 1) + tau(ia, 2)
1158 END DO
1159 rho_set%owns%tau = .true.
1160 rho_set%has%tau = .true.
1161 ELSE
1162 DO ia = 1, na
1163 rho_set%tau(ia, 1, 1) = tau(ia, 1)
1164 END DO
1165 rho_set%owns%tau = .true.
1166 rho_set%has%tau = .true.
1167 END IF
1168 END IF
1169 IF (needs%tau_spin) THEN
1170 cpassert(nspins == 2)
1171 DO ia = 1, na
1172 rho_set%tau_a(ia, 1, 1) = tau(ia, 1)
1173 rho_set%tau_b(ia, 1, 1) = tau(ia, 2)
1174 END DO
1175 rho_set%owns%tau_spin = .true.
1176 rho_set%has%tau_spin = .true.
1177 END IF
1178
1179 ! Laplace
1180 IF (needs%laplace_rho) THEN
1181 IF (nspins == 2) THEN
1182 DO ia = 1, na
1183 rho_set%laplace_rho(ia, 1, 1) = lap(ia, 1) + lap(ia, 2)
1184 END DO
1185 rho_set%owns%laplace_rho = .true.
1186 rho_set%has%laplace_rho = .true.
1187 ELSE
1188 DO ia = 1, na
1189 rho_set%laplace_rho(ia, 1, 1) = lap(ia, 1)
1190 END DO
1191 rho_set%owns%laplace_rho = .true.
1192 rho_set%has%laplace_rho = .true.
1193 END IF
1194 END IF
1195 IF (needs%laplace_rho_spin) THEN
1196 cpassert(nspins == 2)
1197 DO ia = 1, na
1198 rho_set%laplace_rhoa(ia, 1, 1) = lap(ia, 1)
1199 rho_set%laplace_rhob(ia, 1, 1) = lap(ia, 2)
1200 END DO
1201 rho_set%owns%laplace_rho_spin = .true.
1202 rho_set%has%laplace_rho_spin = .true.
1203 END IF
1204
1205 END SUBROUTINE fill_rho_set
1206
1207! **************************************************************************************************
1208!> \brief Calculate PADE exchange-correlation (xc) energy functional and potential
1209!> in spin-restricted case.
1210!> \param rho electron density
1211!> \param exc XC energy functional
1212!> \param vxc XC potential
1213!> \par History
1214!> * 12.2008 created [Juerg Hutter]
1215! **************************************************************************************************
1216 PURE SUBROUTINE lda_pade(rho, exc, vxc)
1217
1218 REAL(dp), DIMENSION(:), INTENT(IN) :: rho
1219 REAL(dp), DIMENSION(:), INTENT(OUT) :: exc, vxc
1220
1221 REAL(kind=dp), PARAMETER :: a0 = 0.4581652932831429e+0_dp, a1 = 0.2217058676663745e+1_dp, &
1222 a2 = 0.7405551735357053e+0_dp, a3 = 0.1968227878617998e-1_dp, &
1223 b1 = 1.0000000000000000e+0_dp, b2 = 0.4504130959426697e+1_dp, &
1224 b3 = 0.1110667363742916e+1_dp, b4 = 0.2359291751427506e-1_dp, f13 = (1.0_dp/3.0_dp), &
1225 rsfac = 0.6203504908994000166680065_dp
1226
1227 INTEGER :: i, n
1228 REAL(kind=dp) :: depade, dpv, dq, epade, p, q, rs
1229
1230 n = SIZE(rho)
1231 exc(1:n) = 0._dp
1232 vxc(1:n) = 0._dp
1233
1234 DO i = 1, n
1235 IF (rho(i) > 1.e-20_dp) THEN
1236 rs = rsfac*rho(i)**(-f13)
1237 p = a0 + (a1 + (a2 + a3*rs)*rs)*rs
1238 q = (b1 + (b2 + (b3 + b4*rs)*rs)*rs)*rs
1239 epade = -p/q
1240
1241 dpv = a1 + (2.0_dp*a2 + 3.0_dp*a3*rs)*rs
1242 dq = b1 + (2.0_dp*b2 + (3.0_dp*b3 + 4.0_dp*b4*rs)*rs)*rs
1243 depade = f13*rs*(dpv*q - p*dq)/(q*q)
1244
1245 exc(i) = epade*rho(i)
1246 vxc(i) = epade + depade
1247 END IF
1248 END DO
1249
1250 END SUBROUTINE lda_pade
1251
1252! **************************************************************************************************
1253!> \brief Calculate PADE exchange-correlation (xc) energy functional and potential
1254!> in spin-polarized case.
1255!> \param rhoa spin-alpha electron density
1256!> \param rhob spin-beta electron density
1257!> \param exc XC energy functional
1258!> \param vxca spin-alpha XC potential
1259!> \param vxcb spin-beta XC potential
1260!> \par History
1261!> * 02.2010 created [Juerg Hutter]
1262! **************************************************************************************************
1263 PURE SUBROUTINE lsd_pade(rhoa, rhob, exc, vxca, vxcb)
1264
1265 REAL(dp), DIMENSION(:), INTENT(IN) :: rhoa, rhob
1266 REAL(dp), DIMENSION(:), INTENT(OUT) :: exc, vxca, vxcb
1267
1268 REAL(kind=dp), PARAMETER :: a0 = 0.4581652932831429e+0_dp, a1 = 0.2217058676663745e+1_dp, &
1269 a2 = 0.7405551735357053e+0_dp, a3 = 0.1968227878617998e-1_dp, &
1270 b1 = 1.0000000000000000e+0_dp, b2 = 0.4504130959426697e+1_dp, &
1271 b3 = 0.1110667363742916e+1_dp, b4 = 0.2359291751427506e-1_dp, &
1272 da0 = 0.119086804055547e+0_dp, da1 = 0.6157402568883345e+0_dp, &
1273 da2 = 0.1574201515892867e+0_dp, da3 = 0.3532336663397157e-2_dp, &
1274 db1 = 0.0000000000000000e+0_dp, db2 = 0.2673612973836267e+0_dp, &
1275 db3 = 0.2052004607777787e+0_dp, db4 = 0.4200005045691381e-2_dp, f13 = (1.0_dp/3.0_dp), &
1276 f43 = (4.0_dp/3.0_dp)
1277 REAL(kind=dp), PARAMETER :: fxfac = 1.923661050931536319759455_dp, &
1278 rsfac = 0.6203504908994000166680065_dp
1279
1280 INTEGER :: i, n
1281 REAL(kind=dp) :: dc, dpv, dq, dr, dx, fa0, fa1, fa2, fa3, &
1282 fb1, fb2, fb3, fb4, fx1, fx2, p, q, &
1283 rhoab, rs, x, xp, xq
1284
1285! 1/(2^(4/3) - 2)
1286
1287 n = SIZE(rhoa)
1288 exc(1:n) = 0._dp
1289 vxca(1:n) = 0._dp
1290 vxcb(1:n) = 0._dp
1291
1292 DO i = 1, n
1293 rhoab = rhoa(i) + rhob(i)
1294 IF (rhoab > 1.e-20_dp) THEN
1295 rs = rsfac*rhoab**(-f13)
1296
1297 x = (rhoa(i) - rhob(i))/rhoab
1298 IF (x < -1.0_dp) THEN
1299 fx1 = 1.0_dp
1300 fx2 = -f43*fxfac*2.0_dp**f13
1301 ELSE IF (x > 1.0_dp) THEN
1302 fx1 = 1.0_dp
1303 fx2 = f43*fxfac*2.0_dp**f13
1304 ELSE
1305 fx1 = ((1.0_dp + x)**f43 + (1.0_dp - x)**f43 - 2.0_dp)*fxfac
1306 fx2 = ((1.0_dp + x)**f13 - (1.0_dp - x)**f13)*fxfac*f43
1307 END IF
1308
1309 fa0 = a0 + fx1*da0
1310 fa1 = a1 + fx1*da1
1311 fa2 = a2 + fx1*da2
1312 fa3 = a3 + fx1*da3
1313 fb1 = b1 + fx1*db1
1314 fb2 = b2 + fx1*db2
1315 fb3 = b3 + fx1*db3
1316 fb4 = b4 + fx1*db4
1317
1318 p = fa0 + (fa1 + (fa2 + fa3*rs)*rs)*rs
1319 q = (fb1 + (fb2 + (fb3 + fb4*rs)*rs)*rs)*rs
1320 dpv = fa1 + (2.0_dp*fa2 + 3.0_dp*fa3*rs)*rs
1321 dq = fb1 + (2.0_dp*fb2 + (3.0_dp*fb3 + &
1322 4.0_dp*fb4*rs)*rs)*rs
1323 xp = da0 + (da1 + (da2 + da3*rs)*rs)*rs
1324 xq = (db1 + (db2 + (db3 + db4*rs)*rs)*rs)*rs
1325
1326 dr = (dpv*q - p*dq)/(q*q)
1327 dx = 2.0_dp*(xp*q - p*xq)/(q*q)*fx2/rhoab
1328 dc = f13*rs*dr - p/q
1329
1330 exc(i) = -p/q*rhoab
1331 vxca(i) = dc - dx*rhob(i)
1332 vxcb(i) = dc + dx*rhoa(i)
1333 END IF
1334 END DO
1335
1336 END SUBROUTINE lsd_pade
1337
1338END MODULE atom_xc
atom_utils::integrate_grid
Definition atom_utils.F:83
atom_types
Define the atom type and its sub types.
Definition atom_types.F:15
atom_utils
Some basic routines for atomic calculations.
Definition atom_utils.F:15
atom_utils::atom_core_density
subroutine, public atom_core_density(corden, potential, typ, rr)
...
Definition atom_utils.F:694
atom_utils::atom_density
subroutine, public atom_density(density, pmat, basis, maxl, typ, rr)
Map the electron density on an atomic radial grid.
Definition atom_utils.F:367
atom_utils::numpot_matrix
subroutine, public numpot_matrix(imat, cpot, basis, derivatives)
Calculate a potential matrix by integrating the potential on an atomic radial grid.
Definition atom_utils.F:1840
atom_utils::coulomb_potential_numeric
subroutine, public coulomb_potential_numeric(cpot, density, grid)
Numerically compute the Coulomb potential on an atomic radial grid.
Definition atom_utils.F:1077
atom_xc
routines that build the integrals of the Vxc potential calculated for the atomic code
Definition atom_xc.F:12
atom_xc::atom_dpot_lda
subroutine, public atom_dpot_lda(basis0, pmat0, basis1, pmat1, maxl, functional, dfexc, linxpar)
Estimate the ADMM exchange energy correction term using a model exchange functional.
Definition atom_xc.F:979
atom_xc::atom_vxc_lsd
subroutine, public atom_vxc_lsd(basis, pmata, pmatb, maxl, xc_section, fexc, xcmata, xcmatb)
Alternative subroutine that calculates atomic exchange-correlation potential in spin-polarized case.
Definition atom_xc.F:800
atom_xc::calculate_atom_vxc_lda
subroutine, public calculate_atom_vxc_lda(xcmat, atom, xc_section)
Calculate atomic exchange-correlation potential in spin-restricted case.
Definition atom_xc.F:199
atom_xc::calculate_atom_ext_vxc
subroutine, public calculate_atom_ext_vxc(vxc, atom, lprint, xcmat)
ZMP subroutine for the scf external density from the external v_xc.
Definition atom_xc.F:148
atom_xc::calculate_atom_zmp
subroutine, public calculate_atom_zmp(ext_density, atom, lprint, xcmat)
ZMP subroutine for the calculation of the effective constraint potential in the atomic code.
Definition atom_xc.F:78
atom_xc::calculate_atom_vxc_lsd
subroutine, public calculate_atom_vxc_lsd(xcmata, xcmatb, atom, xc_section)
Calculate atomic exchange-correlation potential in spin-polarized case.
Definition atom_xc.F:404
atom_xc::atom_vxc_lda
subroutine, public atom_vxc_lda(basis, pmat, maxl, xc_section, fexc, xcmat)
Alternative subroutine that calculates atomic exchange-correlation potential in spin-restricted case.
Definition atom_xc.F:649
atom
Definition atom.F:9
cp_files
Utility routines to open and close files. Tracking of preconnections.
Definition cp_files.F:16
cp_files::open_file
subroutine, public open_file(file_name, file_status, file_form, file_action, file_position, file_pad, unit_number, debug, skip_get_unit_number, file_access)
Opens the requested file using a free unit number.
Definition cp_files.F:308
cp_files::get_unit_number
integer function, public get_unit_number(file_name)
Returns the first logical unit that is not preconnected.
Definition cp_files.F:237
cp_files::close_file
subroutine, public close_file(unit_number, file_status, keep_preconnection)
Close an open file given by its logical unit number. Optionally, keep the file and unit preconnected.
Definition cp_files.F:119
input_constants
collects all constants needed in input so that they can be used without circular dependencies
Definition input_constants.F:17
input_constants::sgp_pseudo
integer, parameter, public sgp_pseudo
Definition input_constants.F:174
input_constants::gth_pseudo
integer, parameter, public gth_pseudo
Definition input_constants.F:174
input_constants::upf_pseudo
integer, parameter, public upf_pseudo
Definition input_constants.F:174
input_constants::xc_none
integer, parameter, public xc_none
Definition input_constants.F:551
input_section_types
objects that represent the structure of input sections and the data contained in an input section
Definition input_section_types.F:15
input_section_types::section_vals_get_subs_vals
recursive type(section_vals_type) function, pointer, public section_vals_get_subs_vals(section_vals, subsection_name, i_rep_section, can_return_null)
returns the values of the requested subsection
Definition input_section_types.F:724
input_section_types::section_vals_val_get
subroutine, public section_vals_val_get(section_vals, keyword_name, i_rep_section, i_rep_val, n_rep_val, val, l_val, i_val, r_val, c_val, l_vals, i_vals, r_vals, c_vals, explicit)
returns the requested value
Definition input_section_types.F:1040
kinds
Defines the basic variable types.
Definition kinds.F:23
kinds::dp
integer, parameter, public dp
Definition kinds.F:34
mathconstants
Definition of mathematical constants and functions.
Definition mathconstants.F:16
mathconstants::pi
real(kind=dp), parameter, public pi
Definition mathconstants.F:110
mathconstants::fourpi
real(kind=dp), parameter, public fourpi
Definition mathconstants.F:113
xc_atom
Definition xc_atom.F:9
xc_atom::xc_rho_set_atom_update
subroutine, public xc_rho_set_atom_update(rho_set, needs, nspins, bo)
...
Definition xc_atom.F:432
xc_derivative_desc
Module with functions to handle derivative descriptors. derivative description are strings have the f...
Definition xc_derivative_desc.F:26
xc_derivative_desc::deriv_norm_drho
integer, parameter, public deriv_norm_drho
Definition xc_derivative_desc.F:35
xc_derivative_desc::deriv_norm_drhoa
integer, parameter, public deriv_norm_drhoa
Definition xc_derivative_desc.F:35
xc_derivative_desc::deriv_rhob
integer, parameter, public deriv_rhob
Definition xc_derivative_desc.F:35
xc_derivative_desc::deriv_rhoa
integer, parameter, public deriv_rhoa
Definition xc_derivative_desc.F:35
xc_derivative_desc::deriv_tau
integer, parameter, public deriv_tau
Definition xc_derivative_desc.F:35
xc_derivative_desc::deriv_tau_b
integer, parameter, public deriv_tau_b
Definition xc_derivative_desc.F:35
xc_derivative_desc::deriv_tau_a
integer, parameter, public deriv_tau_a
Definition xc_derivative_desc.F:35
xc_derivative_desc::deriv_rho
integer, parameter, public deriv_rho
Definition xc_derivative_desc.F:35
xc_derivative_desc::deriv_norm_drhob
integer, parameter, public deriv_norm_drhob
Definition xc_derivative_desc.F:35
xc_derivative_set_types
represent a group ofunctional derivatives
Definition xc_derivative_set_types.F:14
xc_derivative_set_types::xc_dset_zero_all
subroutine, public xc_dset_zero_all(deriv_set)
...
Definition xc_derivative_set_types.F:167
xc_derivative_set_types::xc_dset_get_derivative
type(xc_derivative_type) function, pointer, public xc_dset_get_derivative(derivative_set, description, allocate_deriv)
returns the requested xc_derivative
Definition xc_derivative_set_types.F:66
xc_derivative_set_types::xc_dset_release
subroutine, public xc_dset_release(derivative_set)
releases a derivative set
Definition xc_derivative_set_types.F:144
xc_derivative_set_types::xc_dset_create
subroutine, public xc_dset_create(derivative_set, pw_pool, local_bounds)
creates a derivative set object
Definition xc_derivative_set_types.F:112
xc_derivative_types
Provides types for the management of the xc-functionals and their derivatives.
Definition xc_derivative_types.F:12
xc_derivative_types::xc_derivative_get
subroutine, public xc_derivative_get(deriv, split_desc, order, deriv_data, accept_null_data)
returns various information on the given derivative
Definition xc_derivative_types.F:103
xc_derivatives
Definition xc_derivatives.F:9
xc_derivatives::xc_functionals_get_needs
type(xc_rho_cflags_type) function, public xc_functionals_get_needs(functionals, lsd, calc_potential)
...
Definition xc_derivatives.F:600
xc_derivatives::xc_functionals_eval
subroutine, public xc_functionals_eval(functionals, lsd, rho_set, deriv_set, deriv_order)
...
Definition xc_derivatives.F:565
xc_rho_cflags_types
contains the structure
Definition xc_rho_cflags_types.F:14
xc_rho_set_types
contains the structure
Definition xc_rho_set_types.F:14
xc_rho_set_types::xc_rho_set_create
subroutine, public xc_rho_set_create(rho_set, local_bounds, rho_cutoff, drho_cutoff, tau_cutoff)
allocates and does (minimal) initialization of a rho_set
Definition xc_rho_set_types.F:111
xc_rho_set_types::xc_rho_set_release
subroutine, public xc_rho_set_release(rho_set, pw_pool)
releases the given rho_set
Definition xc_rho_set_types.F:129
atom_types::atom_basis_type
Provides all information about a basis set.
Definition atom_types.F:78
atom_types::atom_type
Provides all information about an atomic kind.
Definition atom_types.F:290
atom_types::opmat_type
Operator matrices.
Definition atom_types.F:248
input_section_types::section_vals_type
stores the values of a section
Definition input_section_types.F:126
xc_derivative_set_types::xc_derivative_set_type
A derivative set contains the different derivatives of a xc-functional in form of a linked list.
Definition xc_derivative_set_types.F:49
xc_derivative_types::xc_derivative_type
represent a derivative of a functional
Definition xc_derivative_types.F:32
xc_rho_cflags_types::xc_rho_cflags_type
contains a flag for each component of xc_rho_set, so that you can use it to tell which components you...
Definition xc_rho_cflags_types.F:48
xc_rho_set_types::xc_rho_set_type
represent a density, with all the representation and data needed to perform a functional evaluation
Definition xc_rho_set_types.F:78