(git:e0f0c17)
Loading...
Searching...
No Matches
qs_vxc.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
10!>
11!>
12!> \par History
13!> refactoring 03-2011 [MI]
14!> \author MI
15! **************************************************************************************************
16MODULE qs_vxc
17
18 USE cell_types, ONLY: cell_type
19 USE cp_control_types, ONLY: dft_control_type
20 USE input_constants, ONLY: sic_ad,&
21 sic_eo,&
22 sic_mauri_spz,&
23 sic_mauri_us,&
24 sic_none,&
25 xc_none,&
26 xc_vdw_fun_nonloc
27 USE input_section_types, ONLY: section_vals_type,&
28 section_vals_val_get
29 USE kinds, ONLY: dp
30 USE message_passing, ONLY: mp_para_env_type
31 USE particle_types, ONLY: particle_type
32 USE pw_env_types, ONLY: pw_env_get,&
33 pw_env_type
34 USE pw_grids, ONLY: pw_grid_compare
35 USE pw_methods, ONLY: pw_axpy,&
36 pw_copy,&
37 pw_multiply,&
38 pw_scale,&
39 pw_transfer,&
40 pw_zero
41 USE pw_pool_types, ONLY: pw_pool_type
42 USE pw_types, ONLY: pw_c1d_gs_type,&
43 pw_r3d_rs_type
44 USE qs_dispersion_nonloc, ONLY: calculate_dispersion_nonloc
45 USE qs_dispersion_types, ONLY: qs_dispersion_type
46 USE qs_ks_types, ONLY: get_ks_env,&
47 qs_ks_env_type
48 USE qs_rho_types, ONLY: qs_rho_get,&
49 qs_rho_type
50 USE skala_gpw_functional, ONLY: skala_gpw_eval,&
51 xc_section_uses_native_skala_grid
52 USE virial_types, ONLY: virial_type
53 USE xc, ONLY: calc_xc_density,&
54 xc_exc_calc,&
55 xc_vxc_pw_create
56#include "./base/base_uses.f90"
57
58 IMPLICIT NONE
59
60 PRIVATE
61
62 ! *** Public subroutines ***
63 PUBLIC :: qs_vxc_create, qs_xc_density
64
65 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'qs_vxc'
66
67CONTAINS
68
69! **************************************************************************************************
70!> \brief calculates and allocates the xc potential, already reducing it to
71!> the dependence on rho and the one on tau
72!> \param ks_env to get all the needed things
73!> \param rho_struct density for which v_xc is calculated
74!> \param xc_section ...
75!> \param vxc_rho will contain the v_xc part that depend on rho
76!> (if one of the chosen xc functionals has it it is allocated and you
77!> are responsible for it)
78!> \param vxc_tau will contain the kinetic tau part of v_xc
79!> (if one of the chosen xc functionals has it it is allocated and you
80!> are responsible for it)
81!> \param exc ...
82!> \param just_energy if true calculates just the energy, and does not
83!> allocate v_*_rspace
84!> \param edisp ...
85!> \param dispersion_env ...
86!> \param adiabatic_rescale_factor ...
87!> \param pw_env_external external plane wave environment
88!> \param native_skala_atom_force ...
89!> \par History
90!> - 05.2002 modified to use the mp_allgather function each pe
91!> computes only part of the grid and this is broadcasted to all
92!> instead of summed.
93!> This scales significantly better (e.g. factor 3 on 12 cpus
94!> 32 H2O) [Joost VdV]
95!> - moved to qs_ks_methods [fawzi]
96!> - sic alterations [Joost VandeVondele]
97!> \author Fawzi Mohamed
98! **************************************************************************************************
99 SUBROUTINE qs_vxc_create(ks_env, rho_struct, xc_section, vxc_rho, vxc_tau, exc, &
100 just_energy, edisp, dispersion_env, adiabatic_rescale_factor, &
101 pw_env_external, native_skala_atom_force)
102
103 TYPE(qs_ks_env_type), POINTER :: ks_env
104 TYPE(qs_rho_type), POINTER :: rho_struct
105 TYPE(section_vals_type), POINTER :: xc_section
106 TYPE(pw_r3d_rs_type), DIMENSION(:), POINTER :: vxc_rho, vxc_tau
107 REAL(kind=dp), INTENT(out) :: exc
108 LOGICAL, INTENT(in), OPTIONAL :: just_energy
109 REAL(kind=dp), INTENT(out), OPTIONAL :: edisp
110 TYPE(qs_dispersion_type), OPTIONAL, POINTER :: dispersion_env
111 REAL(kind=dp), INTENT(in), OPTIONAL :: adiabatic_rescale_factor
112 TYPE(pw_env_type), OPTIONAL, POINTER :: pw_env_external
113 REAL(kind=dp), DIMENSION(:, :), INTENT(OUT), &
114 OPTIONAL :: native_skala_atom_force
115
116 CHARACTER(len=*), PARAMETER :: routinen = 'qs_vxc_create'
117
118 INTEGER :: handle, ispin, mspin, myfun, &
119 nelec_spin(2), vdw
120 LOGICAL :: compute_virial, do_adiabatic_rescaling, my_just_energy, native_skala_grid, &
121 rho_g_valid, sic_scaling_b_zero, tau_g_valid, tau_r_valid, uf_grid, vdw_nl
122 REAL(kind=dp) :: exc_m, factor, &
123 my_adiabatic_rescale_factor, &
124 my_scaling, nelec_s_inv
125 REAL(kind=dp), DIMENSION(3, 3) :: virial_xc_tmp
126 TYPE(cell_type), POINTER :: cell
127 TYPE(dft_control_type), POINTER :: dft_control
128 TYPE(mp_para_env_type), POINTER :: para_env
129 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
130 TYPE(pw_c1d_gs_type), DIMENSION(:), POINTER :: rho_g, rho_m_gspace, rho_struct_g, &
131 tau_struct_g
132 TYPE(pw_c1d_gs_type), POINTER :: rho_nlcc_g, rho_nlcc_g_use, &
133 rho_nlcc_g_xc, tmp_g, tmp_g2
134 TYPE(pw_env_type), POINTER :: pw_env
135 TYPE(pw_pool_type), POINTER :: auxbas_pw_pool, vdw_pw_pool, xc_pw_pool
136 TYPE(pw_r3d_rs_type), DIMENSION(:), POINTER :: my_vxc_rho, my_vxc_tau, rho_m_rspace, &
137 rho_r, rho_struct_r, tau, tau_struct_r
138 TYPE(pw_r3d_rs_type), POINTER :: rho_nlcc, rho_nlcc_use, rho_nlcc_xc, &
139 tmp_pw, weights, weights_use, &
140 weights_xc
141 TYPE(virial_type), POINTER :: virial
142
143 CALL timeset(routinen, handle)
144
145 cpassert(.NOT. ASSOCIATED(vxc_rho))
146 cpassert(.NOT. ASSOCIATED(vxc_tau))
147 NULLIFY (dft_control, pw_env, auxbas_pw_pool, xc_pw_pool, vdw_pw_pool, cell, my_vxc_rho, &
148 tmp_pw, tmp_g, tmp_g2, my_vxc_tau, rho_g, rho_r, tau, rho_m_rspace, &
149 rho_m_gspace, rho_nlcc, rho_nlcc_g, rho_nlcc_g_use, rho_nlcc_g_xc, &
150 rho_nlcc_use, rho_nlcc_xc, rho_struct_r, rho_struct_g, tau_struct_g, tau_struct_r, &
151 weights_use, weights_xc, particle_set)
152
153 exc = 0.0_dp
154 my_just_energy = .false.
155 IF (PRESENT(just_energy)) my_just_energy = just_energy
156 my_adiabatic_rescale_factor = 1.0_dp
157 do_adiabatic_rescaling = .false.
158 IF (PRESENT(adiabatic_rescale_factor)) THEN
159 my_adiabatic_rescale_factor = adiabatic_rescale_factor
160 do_adiabatic_rescaling = .true.
161 END IF
162
163 CALL get_ks_env(ks_env, &
164 dft_control=dft_control, &
165 pw_env=pw_env, &
166 cell=cell, &
167 particle_set=particle_set, &
168 xcint_weights=weights, &
169 virial=virial, &
170 rho_nlcc=rho_nlcc, &
171 rho_nlcc_g=rho_nlcc_g)
172 rho_nlcc_use => rho_nlcc
173 rho_nlcc_g_use => rho_nlcc_g
174 weights_use => weights
175
176 CALL qs_rho_get(rho_struct, &
177 tau_r_valid=tau_r_valid, &
178 tau_g_valid=tau_g_valid, &
179 rho_g_valid=rho_g_valid, &
180 rho_r=rho_struct_r, &
181 rho_g=rho_struct_g, &
182 tau_g=tau_struct_g, &
183 tau_r=tau_struct_r)
184
185 compute_virial = virial%pv_calculate .AND. (.NOT. virial%pv_numer)
186 IF (compute_virial) THEN
187 virial%pv_xc = 0.0_dp
188 END IF
189
190 CALL section_vals_val_get(xc_section, "XC_FUNCTIONAL%_SECTION_PARAMETERS_", &
191 i_val=myfun)
192 CALL section_vals_val_get(xc_section, "VDW_POTENTIAL%POTENTIAL_TYPE", &
193 i_val=vdw)
194
195 vdw_nl = (vdw == xc_vdw_fun_nonloc)
196 ! this combination has not been investigated
197 cpassert(.NOT. (do_adiabatic_rescaling .AND. vdw_nl))
198 ! are the necessary inputs available
199 IF (.NOT. (PRESENT(dispersion_env) .AND. PRESENT(edisp))) THEN
200 vdw_nl = .false.
201 END IF
202 IF (PRESENT(edisp)) edisp = 0.0_dp
203 native_skala_grid = xc_section_uses_native_skala_grid(xc_section)
204 IF (native_skala_grid .AND. ASSOCIATED(rho_nlcc_use)) THEN
205 CALL cp_abort(__location__, &
206 "Native SKALA grid evaluation with NLCC pseudopotentials is not implemented. "// &
207 "The frozen core density would need a SKALA-consistent feature definition.")
208 END IF
209
210 IF (myfun /= xc_none .OR. vdw_nl) THEN
211
212 ! test if the real space density is available
213 cpassert(ASSOCIATED(rho_struct))
214 IF (dft_control%nspins /= 1 .AND. dft_control%nspins /= 2) &
215 cpabort("nspins must be 1 or 2")
216 mspin = SIZE(rho_struct_r)
217 IF (dft_control%nspins == 2 .AND. mspin == 1) &
218 cpabort("Spin count mismatch")
219
220 ! there are some options related to SIC here.
221 ! Normal DFT computes E(rho_alpha,rho_beta) (or its variant E(2*rho_alpha) for non-LSD)
222 ! SIC can E(rho_alpha,rho_beta)-b*(E(rho_alpha,rho_beta)-E(rho_beta,rho_beta))
223 ! or compute E(rho_alpha,rho_beta)-b*E(rho_alpha-rho_beta,0)
224
225 ! my_scaling is the scaling needed of the standard E(rho_alpha,rho_beta) term
226 my_scaling = 1.0_dp
227 SELECT CASE (dft_control%sic_method_id)
228 CASE (sic_none)
229 ! all fine
230 CASE (sic_mauri_spz, sic_ad)
231 ! no idea yet what to do here in that case
232 cpassert(.NOT. tau_r_valid)
233 CASE (sic_mauri_us)
234 my_scaling = 1.0_dp - dft_control%sic_scaling_b
235 ! no idea yet what to do here in that case
236 cpassert(.NOT. tau_r_valid)
237 CASE (sic_eo)
238 ! NOTHING TO BE DONE
239 CASE DEFAULT
240 ! this case has not yet been treated here
241 cpabort("NYI")
242 END SELECT
243
244 IF (dft_control%sic_scaling_b == 0.0_dp) THEN
245 sic_scaling_b_zero = .true.
246 ELSE
247 sic_scaling_b_zero = .false.
248 END IF
249
250 IF (PRESENT(pw_env_external)) &
251 pw_env => pw_env_external
252 CALL pw_env_get(pw_env, xc_pw_pool=xc_pw_pool, auxbas_pw_pool=auxbas_pw_pool)
253 uf_grid = .NOT. pw_grid_compare(auxbas_pw_pool%pw_grid, xc_pw_pool%pw_grid)
254
255 IF (.NOT. uf_grid) THEN
256 rho_r => rho_struct_r
257
258 IF (tau_r_valid) THEN
259 tau => tau_struct_r
260 END IF
261
262 ! for gradient corrected functional the density in g space might
263 ! be useful so if we have it, we pass it in
264 IF (rho_g_valid) THEN
265 rho_g => rho_struct_g
266 END IF
267 ELSE
268 cpassert(rho_g_valid)
269 ALLOCATE (rho_r(mspin))
270 ALLOCATE (rho_g(mspin))
271 DO ispin = 1, mspin
272 CALL xc_pw_pool%create_pw(rho_g(ispin))
273 CALL pw_transfer(rho_struct_g(ispin), rho_g(ispin))
274 END DO
275 DO ispin = 1, mspin
276 CALL xc_pw_pool%create_pw(rho_r(ispin))
277 CALL pw_transfer(rho_g(ispin), rho_r(ispin))
278 END DO
279 IF (tau_r_valid) THEN
280 ALLOCATE (tau(mspin))
281 DO ispin = 1, mspin
282 CALL xc_pw_pool%create_pw(tau(ispin))
283 block
284 TYPE(pw_c1d_gs_type) :: tau_g_aux, tau_g_xc
285 CALL xc_pw_pool%create_pw(tau_g_xc)
286 IF (tau_g_valid) THEN
287 CALL pw_transfer(tau_struct_g(ispin), tau_g_xc)
288 ELSE
289 CALL auxbas_pw_pool%create_pw(tau_g_aux)
290 CALL pw_transfer(tau_struct_r(ispin), tau_g_aux)
291 CALL pw_transfer(tau_g_aux, tau_g_xc)
292 CALL auxbas_pw_pool%give_back_pw(tau_g_aux)
293 END IF
294 CALL pw_transfer(tau_g_xc, tau(ispin))
295 CALL xc_pw_pool%give_back_pw(tau_g_xc)
296 END block
297 END DO
298 END IF
299 IF (ASSOCIATED(weights)) THEN
300 ALLOCATE (weights_xc)
301 CALL xc_pw_pool%create_pw(weights_xc)
302 block
303 TYPE(pw_c1d_gs_type) :: weights_g_aux, weights_g_xc
304 CALL auxbas_pw_pool%create_pw(weights_g_aux)
305 CALL xc_pw_pool%create_pw(weights_g_xc)
306 CALL pw_transfer(weights, weights_g_aux)
307 CALL pw_transfer(weights_g_aux, weights_g_xc)
308 CALL pw_transfer(weights_g_xc, weights_xc)
309 CALL xc_pw_pool%give_back_pw(weights_g_xc)
310 CALL auxbas_pw_pool%give_back_pw(weights_g_aux)
311 END block
312 weights_use => weights_xc
313 END IF
314 IF (ASSOCIATED(rho_nlcc)) THEN
315 cpassert(ASSOCIATED(rho_nlcc_g))
316 ALLOCATE (rho_nlcc_g_xc, rho_nlcc_xc)
317 CALL xc_pw_pool%create_pw(rho_nlcc_g_xc)
318 CALL xc_pw_pool%create_pw(rho_nlcc_xc)
319 CALL pw_transfer(rho_nlcc_g, rho_nlcc_g_xc)
320 CALL pw_transfer(rho_nlcc_g_xc, rho_nlcc_xc)
321 rho_nlcc_use => rho_nlcc_xc
322 rho_nlcc_g_use => rho_nlcc_g_xc
323 END IF
324 END IF
325
326 ! add the nlcc densities
327 IF (ASSOCIATED(rho_nlcc_use)) THEN
328 factor = 1.0_dp
329 DO ispin = 1, mspin
330 CALL pw_axpy(rho_nlcc_use, rho_r(ispin), factor)
331 CALL pw_axpy(rho_nlcc_g_use, rho_g(ispin), factor)
332 END DO
333 END IF
334
335 !
336 ! here the rho_r, rho_g, tau is what it should be
337 ! we get back the right my_vxc_rho and my_vxc_tau as required
338 !
339 IF (native_skala_grid) THEN
340 CALL skala_gpw_eval(vxc_rho=my_vxc_rho, vxc_tau=my_vxc_tau, exc=exc, &
341 rho_r=rho_r, rho_g=rho_g, tau=tau, xc_section=xc_section, &
342 weights=weights_use, pw_pool=xc_pw_pool, &
343 particle_set=particle_set, cell=cell, &
344 compute_virial=compute_virial, virial_xc=virial%pv_xc, &
345 just_energy=my_just_energy, atom_force=native_skala_atom_force)
346 ELSE IF (my_just_energy) THEN
347 exc = xc_exc_calc(rho_r=rho_r, tau=tau, &
348 rho_g=rho_g, xc_section=xc_section, &
349 weights=weights_use, pw_pool=xc_pw_pool)
350
351 ELSE
352 CALL xc_vxc_pw_create(vxc_rho=my_vxc_rho, vxc_tau=my_vxc_tau, rho_r=rho_r, &
353 rho_g=rho_g, tau=tau, exc=exc, &
354 xc_section=xc_section, &
355 weights=weights_use, pw_pool=xc_pw_pool, &
356 compute_virial=compute_virial, &
357 virial_xc=virial%pv_xc)
358 END IF
359
360 ! remove the nlcc densities (keep stuff in original state)
361 IF (ASSOCIATED(rho_nlcc_use)) THEN
362 factor = -1.0_dp
363 DO ispin = 1, mspin
364 CALL pw_axpy(rho_nlcc_use, rho_r(ispin), factor)
365 CALL pw_axpy(rho_nlcc_g_use, rho_g(ispin), factor)
366 END DO
367 END IF
368
369 ! calclulate non-local vdW functional
370 ! only if this XC_SECTION has it
371 ! if yes, we use the dispersion_env from ks_env
372 ! this is dangerous, as it assumes a special connection xc_section -> qs_env
373 IF (vdw_nl) THEN
374 CALL get_ks_env(ks_env=ks_env, para_env=para_env)
375 ! no SIC functionals allowed
376 cpassert(dft_control%sic_method_id == sic_none)
377 !
378 CALL pw_env_get(pw_env, vdw_pw_pool=vdw_pw_pool)
379 IF (my_just_energy) THEN
380 CALL calculate_dispersion_nonloc(my_vxc_rho, rho_r, rho_g, edisp, dispersion_env, &
381 my_just_energy, vdw_pw_pool, xc_pw_pool, para_env)
382 ELSE
383 CALL calculate_dispersion_nonloc(my_vxc_rho, rho_r, rho_g, edisp, dispersion_env, &
384 my_just_energy, vdw_pw_pool, xc_pw_pool, para_env, virial=virial)
385 END IF
386 END IF
387
388 !! Apply rescaling to the potential if requested
389 IF (.NOT. my_just_energy) THEN
390 IF (do_adiabatic_rescaling) THEN
391 IF (ASSOCIATED(my_vxc_rho)) THEN
392 DO ispin = 1, SIZE(my_vxc_rho)
393 CALL pw_scale(my_vxc_rho(ispin), my_adiabatic_rescale_factor)
394 END DO
395 END IF
396 END IF
397 END IF
398
399 IF (my_scaling /= 1.0_dp) THEN
400 exc = exc*my_scaling
401 IF (ASSOCIATED(my_vxc_rho)) THEN
402 DO ispin = 1, SIZE(my_vxc_rho)
403 CALL pw_scale(my_vxc_rho(ispin), my_scaling)
404 END DO
405 END IF
406 IF (ASSOCIATED(my_vxc_tau)) THEN
407 DO ispin = 1, SIZE(my_vxc_tau)
408 CALL pw_scale(my_vxc_tau(ispin), my_scaling)
409 END DO
410 END IF
411 END IF
412
413 ! we have pw data for the xc, qs_ks requests coeff structure, here we transfer
414 ! pw -> coeff
415 IF (ASSOCIATED(my_vxc_rho)) THEN
416 vxc_rho => my_vxc_rho
417 NULLIFY (my_vxc_rho)
418 END IF
419 IF (ASSOCIATED(my_vxc_tau)) THEN
420 vxc_tau => my_vxc_tau
421 NULLIFY (my_vxc_tau)
422 END IF
423 IF (uf_grid) THEN
424 DO ispin = 1, SIZE(rho_r)
425 CALL xc_pw_pool%give_back_pw(rho_r(ispin))
426 END DO
427 DEALLOCATE (rho_r)
428 IF (ASSOCIATED(rho_g)) THEN
429 DO ispin = 1, SIZE(rho_g)
430 CALL xc_pw_pool%give_back_pw(rho_g(ispin))
431 END DO
432 DEALLOCATE (rho_g)
433 END IF
434 END IF
435
436 ! compute again the xc but now for Exc(m,o) and the opposite sign
437 IF (dft_control%sic_method_id == sic_mauri_spz .AND. .NOT. sic_scaling_b_zero) THEN
438 ALLOCATE (rho_m_rspace(2), rho_m_gspace(2))
439 CALL xc_pw_pool%create_pw(rho_m_gspace(1))
440 CALL xc_pw_pool%create_pw(rho_m_rspace(1))
441 CALL pw_copy(rho_struct_r(1), rho_m_rspace(1))
442 CALL pw_axpy(rho_struct_r(2), rho_m_rspace(1), alpha=-1._dp)
443 CALL pw_copy(rho_struct_g(1), rho_m_gspace(1))
444 CALL pw_axpy(rho_struct_g(2), rho_m_gspace(1), alpha=-1._dp)
445 ! bit sad, these will be just zero...
446 CALL xc_pw_pool%create_pw(rho_m_gspace(2))
447 CALL xc_pw_pool%create_pw(rho_m_rspace(2))
448 CALL pw_zero(rho_m_rspace(2))
449 CALL pw_zero(rho_m_gspace(2))
450
451 IF (my_just_energy) THEN
452 exc_m = xc_exc_calc(rho_r=rho_m_rspace, tau=tau, &
453 rho_g=rho_m_gspace, xc_section=xc_section, &
454 weights=weights_use, pw_pool=xc_pw_pool)
455 ELSE
456 ! virial untested
457 cpassert(.NOT. compute_virial)
458 CALL xc_vxc_pw_create(vxc_rho=my_vxc_rho, vxc_tau=my_vxc_tau, rho_r=rho_m_rspace, &
459 rho_g=rho_m_gspace, tau=tau, exc=exc_m, &
460 xc_section=xc_section, &
461 weights=weights_use, pw_pool=xc_pw_pool, &
462 compute_virial=.false., &
463 virial_xc=virial_xc_tmp)
464 END IF
465
466 exc = exc - dft_control%sic_scaling_b*exc_m
467
468 ! and take care of the potential only vxc_rho is taken into account
469 IF (.NOT. my_just_energy) THEN
470 CALL pw_axpy(my_vxc_rho(1), vxc_rho(1), -dft_control%sic_scaling_b)
471 CALL pw_axpy(my_vxc_rho(1), vxc_rho(2), dft_control%sic_scaling_b)
472 CALL my_vxc_rho(1)%release()
473 CALL my_vxc_rho(2)%release()
474 DEALLOCATE (my_vxc_rho)
475 END IF
476
477 DO ispin = 1, 2
478 CALL xc_pw_pool%give_back_pw(rho_m_rspace(ispin))
479 CALL xc_pw_pool%give_back_pw(rho_m_gspace(ispin))
480 END DO
481 DEALLOCATE (rho_m_rspace)
482 DEALLOCATE (rho_m_gspace)
483
484 END IF
485
486 ! now we have - sum_s N_s * Exc(rho_s/N_s,0)
487 IF (dft_control%sic_method_id == sic_ad .AND. .NOT. sic_scaling_b_zero) THEN
488
489 ! find out how many elecs we have
490 CALL get_ks_env(ks_env, nelectron_spin=nelec_spin)
491
492 ALLOCATE (rho_m_rspace(2), rho_m_gspace(2))
493 DO ispin = 1, 2
494 CALL xc_pw_pool%create_pw(rho_m_gspace(ispin))
495 CALL xc_pw_pool%create_pw(rho_m_rspace(ispin))
496 END DO
497
498 DO ispin = 1, 2
499 IF (nelec_spin(ispin) > 0.0_dp) THEN
500 nelec_s_inv = 1.0_dp/nelec_spin(ispin)
501 ELSE
502 ! does it matter if there are no electrons with this spin (H) ?
503 nelec_s_inv = 0.0_dp
504 END IF
505 CALL pw_copy(rho_struct_r(ispin), rho_m_rspace(1))
506 CALL pw_copy(rho_struct_g(ispin), rho_m_gspace(1))
507 CALL pw_scale(rho_m_rspace(1), nelec_s_inv)
508 CALL pw_scale(rho_m_gspace(1), nelec_s_inv)
509 CALL pw_zero(rho_m_rspace(2))
510 CALL pw_zero(rho_m_gspace(2))
511
512 IF (my_just_energy) THEN
513 exc_m = xc_exc_calc(rho_r=rho_m_rspace, tau=tau, &
514 rho_g=rho_m_gspace, xc_section=xc_section, &
515 weights=weights_use, pw_pool=xc_pw_pool)
516 ELSE
517 ! virial untested
518 cpassert(.NOT. compute_virial)
519 CALL xc_vxc_pw_create(vxc_rho=my_vxc_rho, vxc_tau=my_vxc_tau, rho_r=rho_m_rspace, &
520 rho_g=rho_m_gspace, tau=tau, exc=exc_m, &
521 xc_section=xc_section, &
522 weights=weights_use, pw_pool=xc_pw_pool, &
523 compute_virial=.false., &
524 virial_xc=virial_xc_tmp)
525 END IF
526
527 exc = exc - dft_control%sic_scaling_b*nelec_spin(ispin)*exc_m
528
529 ! and take care of the potential only vxc_rho is taken into account
530 IF (.NOT. my_just_energy) THEN
531 CALL pw_axpy(my_vxc_rho(1), vxc_rho(ispin), -dft_control%sic_scaling_b)
532 CALL my_vxc_rho(1)%release()
533 CALL my_vxc_rho(2)%release()
534 DEALLOCATE (my_vxc_rho)
535 END IF
536 END DO
537
538 DO ispin = 1, 2
539 CALL xc_pw_pool%give_back_pw(rho_m_rspace(ispin))
540 CALL xc_pw_pool%give_back_pw(rho_m_gspace(ispin))
541 END DO
542 DEALLOCATE (rho_m_rspace)
543 DEALLOCATE (rho_m_gspace)
544
545 END IF
546
547 ! compute again the xc but now for Exc(n_down,n_down)
548 IF (dft_control%sic_method_id == sic_mauri_us .AND. .NOT. sic_scaling_b_zero) THEN
549 ALLOCATE (rho_r(2))
550 rho_r(1) = rho_struct_r(2)
551 rho_r(2) = rho_struct_r(2)
552 IF (rho_g_valid) THEN
553 ALLOCATE (rho_g(2))
554 rho_g(1) = rho_struct_g(2)
555 rho_g(2) = rho_struct_g(2)
556 END IF
557
558 IF (my_just_energy) THEN
559 exc_m = xc_exc_calc(rho_r=rho_r, tau=tau, &
560 rho_g=rho_g, xc_section=xc_section, &
561 weights=weights_use, pw_pool=xc_pw_pool)
562 ELSE
563 ! virial untested
564 cpassert(.NOT. compute_virial)
565 CALL xc_vxc_pw_create(vxc_rho=my_vxc_rho, vxc_tau=my_vxc_tau, rho_r=rho_r, &
566 rho_g=rho_g, tau=tau, exc=exc_m, &
567 xc_section=xc_section, &
568 weights=weights_use, pw_pool=xc_pw_pool, &
569 compute_virial=.false., &
570 virial_xc=virial_xc_tmp)
571 END IF
572
573 exc = exc + dft_control%sic_scaling_b*exc_m
574
575 ! and take care of the potential
576 IF (.NOT. my_just_energy) THEN
577 ! both go to minority spin
578 CALL pw_axpy(my_vxc_rho(1), vxc_rho(2), 2.0_dp*dft_control%sic_scaling_b)
579 CALL my_vxc_rho(1)%release()
580 CALL my_vxc_rho(2)%release()
581 DEALLOCATE (my_vxc_rho)
582 END IF
583 DEALLOCATE (rho_r, rho_g)
584
585 END IF
586
587 !
588 ! cleanups
589 !
590 IF (uf_grid .AND. (ASSOCIATED(vxc_rho) .OR. ASSOCIATED(vxc_tau))) THEN
591 block
592 TYPE(pw_r3d_rs_type) :: tmp_pw
593 TYPE(pw_c1d_gs_type) :: tmp_g, tmp_g2
594 CALL xc_pw_pool%create_pw(tmp_g)
595 CALL auxbas_pw_pool%create_pw(tmp_g2)
596 IF (ASSOCIATED(vxc_rho)) THEN
597 DO ispin = 1, SIZE(vxc_rho)
598 CALL auxbas_pw_pool%create_pw(tmp_pw)
599 CALL pw_transfer(vxc_rho(ispin), tmp_g)
600 CALL pw_transfer(tmp_g, tmp_g2)
601 CALL pw_transfer(tmp_g2, tmp_pw)
602 CALL xc_pw_pool%give_back_pw(vxc_rho(ispin))
603 vxc_rho(ispin) = tmp_pw
604 END DO
605 END IF
606 IF (ASSOCIATED(vxc_tau)) THEN
607 DO ispin = 1, SIZE(vxc_tau)
608 CALL auxbas_pw_pool%create_pw(tmp_pw)
609 CALL pw_transfer(vxc_tau(ispin), tmp_g)
610 CALL pw_transfer(tmp_g, tmp_g2)
611 CALL pw_transfer(tmp_g2, tmp_pw)
612 CALL xc_pw_pool%give_back_pw(vxc_tau(ispin))
613 vxc_tau(ispin) = tmp_pw
614 END DO
615 END IF
616 CALL auxbas_pw_pool%give_back_pw(tmp_g2)
617 CALL xc_pw_pool%give_back_pw(tmp_g)
618 END block
619 END IF
620 IF (ASSOCIATED(tau) .AND. uf_grid) THEN
621 DO ispin = 1, SIZE(tau)
622 CALL xc_pw_pool%give_back_pw(tau(ispin))
623 END DO
624 DEALLOCATE (tau)
625 END IF
626 IF (ASSOCIATED(weights_xc)) THEN
627 CALL xc_pw_pool%give_back_pw(weights_xc)
628 DEALLOCATE (weights_xc)
629 END IF
630 IF (ASSOCIATED(rho_nlcc_xc)) THEN
631 CALL xc_pw_pool%give_back_pw(rho_nlcc_xc)
632 DEALLOCATE (rho_nlcc_xc)
633 END IF
634 IF (ASSOCIATED(rho_nlcc_g_xc)) THEN
635 CALL xc_pw_pool%give_back_pw(rho_nlcc_g_xc)
636 DEALLOCATE (rho_nlcc_g_xc)
637 END IF
638
639 END IF
640
641 CALL timestop(handle)
642
643 END SUBROUTINE qs_vxc_create
644
645! **************************************************************************************************
646!> \brief calculates the XC density: E_xc(r) - V_xc(r)*rho(r) or E_xc(r)/rho(r)
647!> \param ks_env to get all the needed things
648!> \param rho_struct density
649!> \param xc_section ...
650!> \param dispersion_env ...
651!> \param xc_ener will contain the xc energy density E_xc(r) - V_xc(r)*rho(r)
652!> \param xc_den will contain the xc energy density E_xc(r)/rho(r)
653!> \param exc will contain the xc energy density E_xc(r)
654!> \param vxc ...
655!> \param vtau ...
656!> \author JGH
657! **************************************************************************************************
658 SUBROUTINE qs_xc_density(ks_env, rho_struct, xc_section, dispersion_env, &
659 xc_ener, xc_den, exc, vxc, vtau)
660
661 TYPE(qs_ks_env_type), POINTER :: ks_env
662 TYPE(qs_rho_type), POINTER :: rho_struct
663 TYPE(section_vals_type), POINTER :: xc_section
664 TYPE(qs_dispersion_type), OPTIONAL, POINTER :: dispersion_env
665 TYPE(pw_r3d_rs_type), INTENT(INOUT), OPTIONAL :: xc_ener, xc_den
666 TYPE(pw_r3d_rs_type), OPTIONAL :: exc
667 TYPE(pw_r3d_rs_type), DIMENSION(:), OPTIONAL :: vxc, vtau
668
669 CHARACTER(len=*), PARAMETER :: routinen = 'qs_xc_density'
670
671 INTEGER :: handle, ispin, mspin, myfun, nspins, vdw
672 LOGICAL :: rho_g_valid, tau_g_valid, tau_r_valid, &
673 uf_grid, vdw_nl
674 REAL(kind=dp) :: edisp, excint, factor, rho_cutoff
675 REAL(kind=dp), DIMENSION(3, 3) :: vdum
676 TYPE(cell_type), POINTER :: cell
677 TYPE(dft_control_type), POINTER :: dft_control
678 TYPE(mp_para_env_type), POINTER :: para_env
679 TYPE(pw_c1d_gs_type), DIMENSION(:), POINTER :: rho_g, rho_struct_g, tau_g, tau_struct_g
680 TYPE(pw_c1d_gs_type), POINTER :: rho_nlcc_g, rho_nlcc_g_use, rho_nlcc_g_xc
681 TYPE(pw_env_type), POINTER :: pw_env
682 TYPE(pw_pool_type), POINTER :: auxbas_pw_pool, vdw_pw_pool, xc_pw_pool
683 TYPE(pw_r3d_rs_type) :: exc_r
684 TYPE(pw_r3d_rs_type), DIMENSION(:), POINTER :: rho_r, rho_struct_r, tau_r, &
685 tau_struct_r, vxc_rho, vxc_tau
686 TYPE(pw_r3d_rs_type), POINTER :: rho_nlcc, rho_nlcc_use, rho_nlcc_xc, &
687 weights, weights_use, weights_xc
688
689 CALL timeset(routinen, handle)
690
691 NULLIFY (dft_control, pw_env, auxbas_pw_pool, xc_pw_pool, vdw_pw_pool, cell, &
692 rho_g, rho_struct_g, tau_g, tau_struct_g, rho_nlcc, rho_nlcc_g, &
693 rho_nlcc_g_use, rho_nlcc_g_xc, rho_nlcc_use, rho_nlcc_xc, rho_r, &
694 rho_struct_r, tau_r, tau_struct_r, vxc_rho, vxc_tau, weights, &
695 weights_use, weights_xc)
696
697 CALL get_ks_env(ks_env, &
698 dft_control=dft_control, &
699 pw_env=pw_env, &
700 cell=cell, &
701 xcint_weights=weights, &
702 rho_nlcc=rho_nlcc, &
703 rho_nlcc_g=rho_nlcc_g)
704
705 CALL qs_rho_get(rho_struct, &
706 tau_r_valid=tau_r_valid, &
707 tau_g_valid=tau_g_valid, &
708 rho_g_valid=rho_g_valid, &
709 rho_r=rho_struct_r, &
710 rho_g=rho_struct_g, &
711 tau_r=tau_struct_r, &
712 tau_g=tau_struct_g)
713 nspins = dft_control%nspins
714 mspin = SIZE(rho_struct_r)
715 rho_r => rho_struct_r
716 rho_g => rho_struct_g
717 tau_r => tau_struct_r
718 tau_g => tau_struct_g
719 rho_nlcc_use => rho_nlcc
720 rho_nlcc_g_use => rho_nlcc_g
721 weights_use => weights
722
723 CALL section_vals_val_get(xc_section, "XC_FUNCTIONAL%_SECTION_PARAMETERS_", i_val=myfun)
724 CALL section_vals_val_get(xc_section, "VDW_POTENTIAL%POTENTIAL_TYPE", i_val=vdw)
725 vdw_nl = (vdw == xc_vdw_fun_nonloc)
726 IF (PRESENT(xc_ener)) THEN
727 IF (tau_r_valid) THEN
728 CALL cp_warn(__location__, "Tau contribution will not be correctly handled")
729 END IF
730 END IF
731 IF (vdw_nl) THEN
732 CALL cp_warn(__location__, "vdW functional contribution will be ignored")
733 END IF
734
735 CALL pw_env_get(pw_env, xc_pw_pool=xc_pw_pool, auxbas_pw_pool=auxbas_pw_pool)
736 uf_grid = .NOT. pw_grid_compare(auxbas_pw_pool%pw_grid, xc_pw_pool%pw_grid)
737
738 IF (PRESENT(xc_ener)) THEN
739 CALL pw_zero(xc_ener)
740 END IF
741 IF (PRESENT(xc_den)) THEN
742 CALL pw_zero(xc_den)
743 END IF
744 IF (PRESENT(exc)) THEN
745 CALL pw_zero(exc)
746 END IF
747 IF (PRESENT(vxc)) THEN
748 DO ispin = 1, nspins
749 CALL pw_zero(vxc(ispin))
750 END DO
751 END IF
752 IF (PRESENT(vtau)) THEN
753 DO ispin = 1, nspins
754 CALL pw_zero(vtau(ispin))
755 END DO
756 END IF
757
758 IF (myfun /= xc_none) THEN
759
760 cpassert(ASSOCIATED(rho_struct))
761 cpassert(dft_control%sic_method_id == sic_none)
762
763 IF (uf_grid) THEN
764 NULLIFY (rho_r, rho_g, tau_r, tau_g)
765 IF (rho_g_valid) THEN
766 CALL create_density_on_pool(xc_pw_pool, rho_struct_g, rho_r, rho_g)
767 ELSEIF (ASSOCIATED(rho_struct_r)) THEN
768 CALL create_density_on_pool_from_r(auxbas_pw_pool, xc_pw_pool, rho_struct_r, rho_r, rho_g)
769 ELSE
770 cpabort("Fine Grid in qs_xc_density requires rho_r or rho_g")
771 END IF
772 IF (tau_r_valid) THEN
773 IF (tau_g_valid) THEN
774 CALL create_density_on_pool(xc_pw_pool, tau_struct_g, tau_r, tau_g)
775 ELSEIF (ASSOCIATED(tau_struct_r)) THEN
776 CALL create_density_on_pool_from_r(auxbas_pw_pool, xc_pw_pool, tau_struct_r, tau_r, tau_g)
777 ELSE
778 cpabort("Fine Grid in qs_xc_density requires tau_r or tau_g")
779 END IF
780 END IF
781 IF (ASSOCIATED(weights)) THEN
782 ALLOCATE (weights_xc)
783 CALL xc_pw_pool%create_pw(weights_xc)
784 CALL transfer_rspace_between_pools(auxbas_pw_pool, xc_pw_pool, weights, weights_xc)
785 weights_use => weights_xc
786 END IF
787 IF (ASSOCIATED(rho_nlcc)) THEN
788 cpassert(ASSOCIATED(rho_nlcc_g))
789 ALLOCATE (rho_nlcc_g_xc, rho_nlcc_xc)
790 CALL xc_pw_pool%create_pw(rho_nlcc_g_xc)
791 CALL xc_pw_pool%create_pw(rho_nlcc_xc)
792 CALL pw_transfer(rho_nlcc_g, rho_nlcc_g_xc)
793 CALL pw_transfer(rho_nlcc_g_xc, rho_nlcc_xc)
794 rho_nlcc_use => rho_nlcc_xc
795 rho_nlcc_g_use => rho_nlcc_g_xc
796 END IF
797 END IF
798
799 ! add the nlcc densities
800 IF (ASSOCIATED(rho_nlcc_use)) THEN
801 factor = 1.0_dp
802 DO ispin = 1, mspin
803 CALL pw_axpy(rho_nlcc_use, rho_r(ispin), factor)
804 CALL pw_axpy(rho_nlcc_g_use, rho_g(ispin), factor)
805 END DO
806 END IF
807 NULLIFY (vxc_rho, vxc_tau)
808 CALL xc_vxc_pw_create(vxc_rho=vxc_rho, vxc_tau=vxc_tau, rho_r=rho_r, &
809 rho_g=rho_g, tau=tau_r, exc=excint, &
810 xc_section=xc_section, &
811 weights=weights_use, pw_pool=xc_pw_pool, &
812 compute_virial=.false., &
813 virial_xc=vdum, &
814 exc_r=exc_r)
815 ! calclulate non-local vdW functional
816 ! only if this XC_SECTION has it
817 ! if yes, we use the dispersion_env from ks_env
818 ! this is dangerous, as it assumes a special connection xc_section -> qs_env
819 IF (vdw_nl) THEN
820 CALL get_ks_env(ks_env=ks_env, para_env=para_env)
821 ! no SIC functionals allowed
822 cpassert(dft_control%sic_method_id == sic_none)
823 !
824 CALL pw_env_get(pw_env, vdw_pw_pool=vdw_pw_pool)
825 CALL calculate_dispersion_nonloc(vxc_rho, rho_r, rho_g, edisp, dispersion_env, &
826 .false., vdw_pw_pool, xc_pw_pool, para_env)
827 END IF
828
829 ! remove the nlcc densities (keep stuff in original state)
830 IF (ASSOCIATED(rho_nlcc_use)) THEN
831 factor = -1.0_dp
832 DO ispin = 1, mspin
833 CALL pw_axpy(rho_nlcc_use, rho_r(ispin), factor)
834 CALL pw_axpy(rho_nlcc_g_use, rho_g(ispin), factor)
835 END DO
836 END IF
837 !
838 IF (PRESENT(xc_den)) THEN
839 rho_cutoff = 1.e-14_dp
840 IF (uf_grid) THEN
841 block
842 TYPE(pw_r3d_rs_type) :: tmp_pw
843 CALL xc_pw_pool%create_pw(tmp_pw)
844 CALL pw_copy(exc_r, tmp_pw)
845 CALL calc_xc_density(tmp_pw, rho_r, rho_cutoff)
846 CALL transfer_rspace_between_pools(xc_pw_pool, auxbas_pw_pool, tmp_pw, xc_den)
847 CALL xc_pw_pool%give_back_pw(tmp_pw)
848 END block
849 ELSE
850 CALL pw_copy(exc_r, xc_den)
851 CALL calc_xc_density(xc_den, rho_r, rho_cutoff)
852 END IF
853 END IF
854 IF (PRESENT(xc_ener)) THEN
855 IF (uf_grid) THEN
856 block
857 TYPE(pw_r3d_rs_type) :: tmp_pw
858 CALL xc_pw_pool%create_pw(tmp_pw)
859 CALL pw_copy(exc_r, tmp_pw)
860 DO ispin = 1, nspins
861 CALL pw_multiply(tmp_pw, vxc_rho(ispin), rho_r(ispin), alpha=-1.0_dp)
862 END DO
863 CALL transfer_rspace_between_pools(xc_pw_pool, auxbas_pw_pool, tmp_pw, xc_ener)
864 CALL xc_pw_pool%give_back_pw(tmp_pw)
865 END block
866 ELSE
867 CALL pw_copy(exc_r, xc_ener)
868 DO ispin = 1, nspins
869 CALL pw_multiply(xc_ener, vxc_rho(ispin), rho_r(ispin), alpha=-1.0_dp)
870 END DO
871 END IF
872 END IF
873 IF (PRESENT(exc)) THEN
874 IF (uf_grid) THEN
875 CALL transfer_rspace_between_pools(xc_pw_pool, auxbas_pw_pool, exc_r, exc)
876 ELSE
877 CALL pw_copy(exc_r, exc)
878 END IF
879 END IF
880 IF (PRESENT(vxc)) THEN
881 DO ispin = 1, nspins
882 IF (uf_grid) THEN
883 CALL transfer_rspace_between_pools(xc_pw_pool, auxbas_pw_pool, vxc_rho(ispin), vxc(ispin))
884 ELSE
885 CALL pw_copy(vxc_rho(ispin), vxc(ispin))
886 END IF
887 END DO
888 END IF
889 IF (PRESENT(vtau) .AND. ASSOCIATED(vxc_tau)) THEN
890 DO ispin = 1, nspins
891 IF (uf_grid) THEN
892 CALL transfer_rspace_between_pools(xc_pw_pool, auxbas_pw_pool, vxc_tau(ispin), vtau(ispin))
893 ELSE
894 CALL pw_copy(vxc_tau(ispin), vtau(ispin))
895 END IF
896 END DO
897 END IF
898 ! remove arrays
899 IF (ASSOCIATED(vxc_rho)) THEN
900 DO ispin = 1, nspins
901 CALL vxc_rho(ispin)%release()
902 END DO
903 DEALLOCATE (vxc_rho)
904 END IF
905 IF (ASSOCIATED(vxc_tau)) THEN
906 DO ispin = 1, nspins
907 CALL vxc_tau(ispin)%release()
908 END DO
909 DEALLOCATE (vxc_tau)
910 END IF
911 CALL exc_r%release()
912 IF (uf_grid) THEN
913 CALL give_back_density_on_pool(xc_pw_pool, rho_r, rho_g)
914 IF (ASSOCIATED(tau_r)) CALL give_back_density_on_pool(xc_pw_pool, tau_r, tau_g)
915 IF (ASSOCIATED(weights_xc)) THEN
916 CALL xc_pw_pool%give_back_pw(weights_xc)
917 DEALLOCATE (weights_xc)
918 END IF
919 IF (ASSOCIATED(rho_nlcc_xc)) THEN
920 CALL xc_pw_pool%give_back_pw(rho_nlcc_xc)
921 DEALLOCATE (rho_nlcc_xc)
922 END IF
923 IF (ASSOCIATED(rho_nlcc_g_xc)) THEN
924 CALL xc_pw_pool%give_back_pw(rho_nlcc_g_xc)
925 DEALLOCATE (rho_nlcc_g_xc)
926 END IF
927 END IF
928 !
929 END IF
930
931 CALL timestop(handle)
932
933 END SUBROUTINE qs_xc_density
934
935! **************************************************************************************************
936!> \brief transfers an r-space PW between two pools and writes into an existing target PW
937!> \param source_pw_pool ...
938!> \param target_pw_pool ...
939!> \param source ...
940!> \param TARGET ...
941! **************************************************************************************************
942 SUBROUTINE transfer_rspace_between_pools(source_pw_pool, target_pw_pool, source, TARGET)
943 TYPE(pw_pool_type), POINTER :: source_pw_pool, target_pw_pool
944 TYPE(pw_r3d_rs_type), INTENT(INOUT) :: source, target
945
946 TYPE(pw_c1d_gs_type) :: source_g, target_g
947
948 cpassert(ASSOCIATED(source_pw_pool))
949 cpassert(ASSOCIATED(target_pw_pool))
950
951 IF (pw_grid_compare(source_pw_pool%pw_grid, target_pw_pool%pw_grid)) THEN
952 CALL pw_copy(source, TARGET)
953 ELSE
954 CALL source_pw_pool%create_pw(source_g)
955 CALL target_pw_pool%create_pw(target_g)
956 CALL pw_transfer(source, source_g)
957 CALL pw_transfer(source_g, target_g)
958 CALL pw_transfer(target_g, TARGET)
959 CALL target_pw_pool%give_back_pw(target_g)
960 CALL source_pw_pool%give_back_pw(source_g)
961 END IF
962
963 END SUBROUTINE transfer_rspace_between_pools
964
965! **************************************************************************************************
966!> \brief transfers a g-space density to a given PW pool and creates its r-space representation
967!> \param pw_pool ...
968!> \param rho_g_in ...
969!> \param rho_r_out ...
970!> \param rho_g_out ...
971! **************************************************************************************************
972 SUBROUTINE create_density_on_pool(pw_pool, rho_g_in, rho_r_out, rho_g_out)
973 TYPE(pw_pool_type), POINTER :: pw_pool
974 TYPE(pw_c1d_gs_type), DIMENSION(:), POINTER :: rho_g_in
975 TYPE(pw_r3d_rs_type), DIMENSION(:), POINTER :: rho_r_out
976 TYPE(pw_c1d_gs_type), DIMENSION(:), POINTER :: rho_g_out
977
978 INTEGER :: ispin, nspins
979
980 cpassert(ASSOCIATED(pw_pool))
981 cpassert(ASSOCIATED(rho_g_in))
982
983 nspins = SIZE(rho_g_in)
984 ALLOCATE (rho_r_out(nspins), rho_g_out(nspins))
985 DO ispin = 1, nspins
986 CALL pw_pool%create_pw(rho_g_out(ispin))
987 CALL pw_pool%create_pw(rho_r_out(ispin))
988 CALL pw_transfer(rho_g_in(ispin), rho_g_out(ispin))
989 CALL pw_transfer(rho_g_out(ispin), rho_r_out(ispin))
990 END DO
991
992 END SUBROUTINE create_density_on_pool
993
994! **************************************************************************************************
995!> \brief transfers an r-space density to a given PW pool and creates its g-space representation
996!> \param source_pw_pool ...
997!> \param target_pw_pool ...
998!> \param rho_r_in ...
999!> \param rho_r_out ...
1000!> \param rho_g_out ...
1001! **************************************************************************************************
1002 SUBROUTINE create_density_on_pool_from_r(source_pw_pool, target_pw_pool, rho_r_in, rho_r_out, rho_g_out)
1003 TYPE(pw_pool_type), POINTER :: source_pw_pool, target_pw_pool
1004 TYPE(pw_r3d_rs_type), DIMENSION(:), POINTER :: rho_r_in, rho_r_out
1005 TYPE(pw_c1d_gs_type), DIMENSION(:), POINTER :: rho_g_out
1006
1007 INTEGER :: ispin, nspins
1008 TYPE(pw_c1d_gs_type) :: rho_g_in
1009
1010 cpassert(ASSOCIATED(source_pw_pool))
1011 cpassert(ASSOCIATED(target_pw_pool))
1012 cpassert(ASSOCIATED(rho_r_in))
1013
1014 nspins = SIZE(rho_r_in)
1015 ALLOCATE (rho_r_out(nspins), rho_g_out(nspins))
1016 DO ispin = 1, nspins
1017 CALL source_pw_pool%create_pw(rho_g_in)
1018 CALL target_pw_pool%create_pw(rho_g_out(ispin))
1019 CALL target_pw_pool%create_pw(rho_r_out(ispin))
1020 CALL pw_transfer(rho_r_in(ispin), rho_g_in)
1021 CALL pw_transfer(rho_g_in, rho_g_out(ispin))
1022 CALL pw_transfer(rho_g_out(ispin), rho_r_out(ispin))
1023 CALL source_pw_pool%give_back_pw(rho_g_in)
1024 END DO
1025
1026 END SUBROUTINE create_density_on_pool_from_r
1027
1028! **************************************************************************************************
1029!> \brief returns temporary density arrays to the given PW pool
1030!> \param pw_pool ...
1031!> \param rho_r ...
1032!> \param rho_g ...
1033! **************************************************************************************************
1034 SUBROUTINE give_back_density_on_pool(pw_pool, rho_r, rho_g)
1035 TYPE(pw_pool_type), POINTER :: pw_pool
1036 TYPE(pw_r3d_rs_type), DIMENSION(:), POINTER :: rho_r
1037 TYPE(pw_c1d_gs_type), DIMENSION(:), POINTER :: rho_g
1038
1039 INTEGER :: ispin
1040
1041 cpassert(ASSOCIATED(pw_pool))
1042
1043 IF (ASSOCIATED(rho_r)) THEN
1044 DO ispin = 1, SIZE(rho_r)
1045 CALL pw_pool%give_back_pw(rho_r(ispin))
1046 END DO
1047 DEALLOCATE (rho_r)
1048 END IF
1049 IF (ASSOCIATED(rho_g)) THEN
1050 DO ispin = 1, SIZE(rho_g)
1051 CALL pw_pool%give_back_pw(rho_g(ispin))
1052 END DO
1053 DEALLOCATE (rho_g)
1054 END IF
1055
1056 END SUBROUTINE give_back_density_on_pool
1057
1058END MODULE qs_vxc
pw_methods::pw_axpy
Definition pw_methods.F:165
pw_methods::pw_copy
Definition pw_methods.F:146
pw_methods::pw_multiply
Definition pw_methods.F:184
pw_methods::pw_scale
Definition pw_methods.F:94
pw_methods::pw_transfer
Definition pw_methods.F:304
pw_methods::pw_zero
Definition pw_methods.F:83
cell_types
Handles all functions related to the CELL.
Definition cell_types.F:15
cp_control_types
Defines control structures, which contain the parameters and the settings for the DFT-based calculati...
Definition cp_control_types.F:12
input_constants
collects all constants needed in input so that they can be used without circular dependencies
Definition input_constants.F:17
input_constants::sic_mauri_spz
integer, parameter, public sic_mauri_spz
Definition input_constants.F:588
input_constants::xc_vdw_fun_nonloc
integer, parameter, public xc_vdw_fun_nonloc
Definition input_constants.F:615
input_constants::sic_eo
integer, parameter, public sic_eo
Definition input_constants.F:588
input_constants::sic_mauri_us
integer, parameter, public sic_mauri_us
Definition input_constants.F:588
input_constants::sic_none
integer, parameter, public sic_none
Definition input_constants.F:588
input_constants::xc_none
integer, parameter, public xc_none
Definition input_constants.F:571
input_constants::sic_ad
integer, parameter, public sic_ad
Definition input_constants.F:588
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_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:1047
kinds
Defines the basic variable types.
Definition kinds.F:23
kinds::dp
integer, parameter, public dp
Definition kinds.F:34
message_passing
Interface to the message passing library MPI.
Definition message_passing.F:23
particle_types
Define the data structure for the particle information.
Definition particle_types.F:19
pw_env_types
container for various plainwaves related things
Definition pw_env_types.F:14
pw_env_types::pw_env_get
subroutine, public pw_env_get(pw_env, pw_pools, cube_info, gridlevel_info, auxbas_pw_pool, auxbas_grid, auxbas_rs_desc, auxbas_rs_grid, rs_descs, rs_grids, xc_pw_pool, vdw_pw_pool, poisson_env, interp_section)
returns the various attributes of the pw env
Definition pw_env_types.F:113
pw_grids
This module defines the grid data type and some basic operations on it.
Definition pw_grids.F:36
pw_grids::pw_grid_compare
logical function, public pw_grid_compare(grida, gridb)
Check if two pw_grids are equal.
Definition pw_grids.F:148
pw_methods
Definition pw_methods.F:25
pw_pool_types
Manages a pool of grids (to be used for example as tmp objects), but can also be used to instantiate ...
Definition pw_pool_types.F:24
pw_types
Definition pw_types.F:24
qs_dispersion_nonloc
Calculation of non local dispersion functionals Some routines adapted from: Copyright (C) 2001-2009 Q...
Definition qs_dispersion_nonloc.F:19
qs_dispersion_nonloc::calculate_dispersion_nonloc
subroutine, public calculate_dispersion_nonloc(vxc_rho, rho_r, rho_g, edispersion, dispersion_env, energy_only, pw_pool, xc_pw_pool, para_env, virial)
Calculates the non-local vdW functional using the method of Soler For spin polarized cases we use E(a...
Definition qs_dispersion_nonloc.F:166
qs_dispersion_types
Definition of disperson types for DFT calculations.
Definition qs_dispersion_types.F:12
qs_ks_types
Definition qs_ks_types.F:15
qs_ks_types::get_ks_env
subroutine, public get_ks_env(ks_env, v_hartree_rspace, s_mstruct_changed, rho_changed, exc_accint, potential_changed, forces_up_to_date, complex_ks, matrix_h, matrix_h_im, matrix_ks, matrix_ks_im, matrix_vxc, kinetic, matrix_s, matrix_s_ri_aux, matrix_w, matrix_p_mp2, matrix_p_mp2_admm, matrix_h_kp, matrix_h_im_kp, matrix_ks_kp, matrix_vxc_kp, kinetic_kp, matrix_s_kp, matrix_w_kp, matrix_s_ri_aux_kp, matrix_ks_im_kp, rho, rho_xc, vppl, xcint_weights, rho_core, rho_nlcc, rho_nlcc_g, vee, neighbor_list_id, sab_orb, sab_all, sac_ae, sac_ppl, sac_lri, sap_ppnl, sap_oce, sab_lrc, sab_se, sab_xtbe, sab_tbe, sab_core, sab_xb, sab_xtb_pp, sab_xtb_nonbond, sab_vdw, sab_scp, sab_almo, sab_kp, sab_kp_nosym, sab_cneo, task_list, task_list_soft, kpoints, do_kpoints, atomic_kind_set, qs_kind_set, cell, cell_ref, use_ref_cell, particle_set, energy, force, local_particles, local_molecules, molecule_kind_set, molecule_set, subsys, cp_subsys, virial, results, atprop, nkind, natom, dft_control, dbcsr_dist, distribution_2d, pw_env, para_env, blacs_env, nelectron_total, nelectron_spin)
...
Definition qs_ks_types.F:341
qs_rho_types
superstucture that hold various representations of the density and keeps track of which ones are vali...
Definition qs_rho_types.F:18
qs_rho_types::qs_rho_get
subroutine, public qs_rho_get(rho_struct, rho_ao, rho_ao_im, rho_ao_kp, rho_ao_im_kp, rho_r, drho_r, rho_g, drho_g, tau_r, tau_g, rho_r_valid, drho_r_valid, rho_g_valid, drho_g_valid, tau_r_valid, tau_g_valid, tot_rho_r, tot_rho_g, rho_r_sccs, soft_valid, complex_rho_ao)
returns info about the density described by this object. If some representation is not available an e...
Definition qs_rho_types.F:229
qs_vxc
Definition qs_vxc.F:16
qs_vxc::qs_vxc_create
subroutine, public qs_vxc_create(ks_env, rho_struct, xc_section, vxc_rho, vxc_tau, exc, just_energy, edisp, dispersion_env, adiabatic_rescale_factor, pw_env_external, native_skala_atom_force)
calculates and allocates the xc potential, already reducing it to the dependence on rho and the one o...
Definition qs_vxc.F:102
qs_vxc::qs_xc_density
subroutine, public qs_xc_density(ks_env, rho_struct, xc_section, dispersion_env, xc_ener, xc_den, exc, vxc, vtau)
calculates the XC density: E_xc(r) - V_xc(r)*rho(r) or E_xc(r)/rho(r)
Definition qs_vxc.F:660
skala_gpw_functional
Experimental CP2K-native GPW real-space-grid path for SKALA TorchScript models.
Definition skala_gpw_functional.F:11
skala_gpw_functional::skala_gpw_eval
subroutine, public skala_gpw_eval(vxc_rho, vxc_tau, exc, rho_r, rho_g, tau, xc_section, weights, pw_pool, particle_set, cell, compute_virial, virial_xc, just_energy, atom_force)
Evaluate SKALA energy and first derivatives on a CP2K GPW grid.
Definition skala_gpw_functional.F:174
skala_gpw_functional::xc_section_uses_native_skala_grid
logical function, public xc_section_uses_native_skala_grid(xc_section)
Return true if the GAUXC subsection requests the CP2K-native GPW grid path.
Definition skala_gpw_functional.F:81
virial_types
Definition virial_types.F:13
xc
Exchange and Correlation functional calculations.
Definition xc.F:17
xc::xc_exc_calc
real(kind=dp) function, public xc_exc_calc(rho_r, rho_g, tau, xc_section, weights, pw_pool)
calculates just the exchange and correlation energy (no vxc)
Definition xc.F:791
xc::xc_vxc_pw_create
subroutine, public xc_vxc_pw_create(vxc_rho, vxc_tau, exc, rho_r, rho_g, tau, xc_section, weights, pw_pool, compute_virial, virial_xc, exc_r)
Exchange and Correlation functional calculations.
Definition xc.F:474
xc::calc_xc_density
subroutine, public calc_xc_density(pot, rho, rho_cutoff)
Definition xc.F:402
cell_types::cell_type
Type defining parameters related to the simulation cell.
Definition cell_types.F:60
cp_control_types::dft_control_type
Definition cp_control_types.F:743
input_section_types::section_vals_type
stores the values of a section
Definition input_section_types.F:127
message_passing::mp_para_env_type
stores all the informations relevant to an mpi environment
Definition message_passing.F:743
particle_types::particle_type
Definition particle_types.F:35
pw_env_types::pw_env_type
contained for different pw related things
Definition pw_env_types.F:70
pw_pool_types::pw_pool_type
Manages a pool of grids (to be used for example as tmp objects), but can also be used to instantiate ...
Definition pw_pool_types.F:83
pw_types::pw_c1d_gs_type
Definition pw_types.F:127
pw_types::pw_r3d_rs_type
Definition pw_types.F:62
qs_dispersion_types::qs_dispersion_type
Definition qs_dispersion_types.F:34
qs_ks_types::qs_ks_env_type
calculation environment to calculate the ks matrix, holds all the needed vars. assumes that the core ...
Definition qs_ks_types.F:137
qs_rho_types::qs_rho_type
keeps the density in various representations, keeping track of which ones are valid.
Definition qs_rho_types.F:62
virial_types::virial_type
Definition virial_types.F:26