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