(git:5e7fe52)
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 cp_dbcsr_api, ONLY: dbcsr_p_type
21 USE cp_log_handling, ONLY: cp_logger_get_default_io_unit
22 USE gapw_gspace_reconstruction, ONLY: calculate_rhotot_elec_gspace
23 USE input_constants, ONLY: sic_ad,&
24 sic_eo,&
25 sic_mauri_spz,&
26 sic_mauri_us,&
27 sic_none,&
28 xc_none,&
29 xc_vdw_fun_nonloc
30 USE input_section_types, ONLY: section_vals_type,&
31 section_vals_val_get
32 USE kinds, ONLY: dp
33 USE message_passing, ONLY: mp_para_env_type
34 USE particle_types, ONLY: particle_type
35 USE pw_env_types, ONLY: pw_env_get,&
36 pw_env_type
37 USE pw_grids, ONLY: get_pw_grid_info,&
38 pw_grid_compare
39 USE pw_methods, ONLY: pw_axpy,&
40 pw_copy,&
41 pw_integral_ab,&
42 pw_integrate_function,&
43 pw_multiply,&
44 pw_scale,&
45 pw_transfer,&
46 pw_zero
47 USE pw_pool_types, ONLY: pw_pool_type
48 USE pw_types, ONLY: pw_c1d_gs_type,&
49 pw_r3d_rs_type
50 USE qs_collocate_density, ONLY: calculate_rho_elec
51 USE qs_dispersion_nonloc, ONLY: calculate_dispersion_nonloc
52 USE qs_dispersion_types, ONLY: qs_dispersion_type
53 USE qs_environment_types, ONLY: qs_environment_type
54 USE qs_ks_types, ONLY: get_ks_env,&
55 qs_ks_env_type
56 USE qs_rho_types, ONLY: qs_rho_get,&
57 qs_rho_type
58 USE skala_gpw_functional, ONLY: get_gauxc_section,&
59 native_skala_gapw_composite_direct_ao,&
60 native_skala_gapw_composite_reference,&
61 skala_gpw_eval,&
62 xc_section_uses_native_skala_evaluator
63 USE virial_types, ONLY: virial_type
64 USE xc, ONLY: calc_xc_density,&
65 xc_exc_calc,&
66 xc_vxc_pw_create
67#include "./base/base_uses.f90"
68
69 IMPLICIT NONE
70
71 PRIVATE
72
73 ! *** Public subroutines ***
74 PUBLIC :: qs_vxc_create, qs_xc_density
75
76 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'qs_vxc'
77
78CONTAINS
79
80! **************************************************************************************************
81!> \brief calculates and allocates the xc potential, already reducing it to
82!> the dependence on rho and the one on tau
83!> \param ks_env to get all the needed things
84!> \param rho_struct density for which v_xc is calculated
85!> \param xc_section ...
86!> \param vxc_rho will contain the v_xc part that depend on rho
87!> (if one of the chosen xc functionals has it it is allocated and you
88!> are responsible for it)
89!> \param vxc_tau will contain the kinetic tau part of v_xc
90!> (if one of the chosen xc functionals has it it is allocated and you
91!> are responsible for it)
92!> \param exc ...
93!> \param just_energy if true calculates just the energy, and does not
94!> allocate v_*_rspace
95!> \param edisp ...
96!> \param dispersion_env ...
97!> \param adiabatic_rescale_factor ...
98!> \param pw_env_external external plane wave environment
99!> \param native_skala_atom_force ...
100!> \param qs_env_external parent Quickstep environment for GAPW reconstruction
101!> \param native_gapw_composite_override ...
102!> \param native_skala_defer_to_atom_composite prepare VXC storage while deferring the
103!> native SKALA evaluation to the atom-centered GAPW composite path
104!> \par History
105!> - 05.2002 modified to use the mp_allgather function each pe
106!> computes only part of the grid and this is broadcasted to all
107!> instead of summed.
108!> This scales significantly better (e.g. factor 3 on 12 cpus
109!> 32 H2O) [Joost VdV]
110!> - moved to qs_ks_methods [fawzi]
111!> - sic alterations [Joost VandeVondele]
112!> \author Fawzi Mohamed
113! **************************************************************************************************
114 SUBROUTINE qs_vxc_create(ks_env, rho_struct, xc_section, vxc_rho, vxc_tau, exc, &
115 just_energy, edisp, dispersion_env, adiabatic_rescale_factor, &
116 pw_env_external, native_skala_atom_force, qs_env_external, &
117 native_gapw_composite_override, native_skala_defer_to_atom_composite)
118
119 TYPE(qs_ks_env_type), POINTER :: ks_env
120 TYPE(qs_rho_type), POINTER :: rho_struct
121 TYPE(section_vals_type), POINTER :: xc_section
122 TYPE(pw_r3d_rs_type), DIMENSION(:), POINTER :: vxc_rho, vxc_tau
123 REAL(kind=dp), INTENT(out) :: exc
124 LOGICAL, INTENT(in), OPTIONAL :: just_energy
125 REAL(kind=dp), INTENT(out), OPTIONAL :: edisp
126 TYPE(qs_dispersion_type), OPTIONAL, POINTER :: dispersion_env
127 REAL(kind=dp), INTENT(in), OPTIONAL :: adiabatic_rescale_factor
128 TYPE(pw_env_type), OPTIONAL, POINTER :: pw_env_external
129 REAL(kind=dp), DIMENSION(:, :), INTENT(OUT), &
130 OPTIONAL :: native_skala_atom_force
131 TYPE(qs_environment_type), OPTIONAL, POINTER :: qs_env_external
132 LOGICAL, INTENT(in), OPTIONAL :: native_gapw_composite_override, &
133 native_skala_defer_to_atom_composite
134
135 CHARACTER(len=*), PARAMETER :: routinen = 'qs_vxc_create'
136
137 INTEGER :: handle, i, ispin, j, k, mspin, myfun, &
138 nelec_spin(2), output_unit, vdw
139 LOGICAL :: compute_virial, defer_native_skala_to_atom_composite, do_adiabatic_rescaling, &
140 my_just_energy, native_gapw_composite_direct_ao, native_gapw_composite_reference, &
141 native_grid_diagnostics, native_skala_grid, rho_g_valid, sic_scaling_b_zero, tau_g_valid, &
142 tau_r_valid, uf_grid, vdw_nl
143 REAL(kind=dp) :: composite_hard_integral, composite_rho_max, composite_rho_min, &
144 composite_rho_r_integral, composite_soft_integral, composite_tau_max, composite_tau_min, &
145 composite_tau_r_integral, delta, direct_rho_integral, direct_tau_integral, exc_m, factor, &
146 my_adiabatic_rescale_factor, my_scaling, nelec_s_inv, q_max, rho_composite_hard, &
147 rho_composite_soft, rho_diff_l2, rho_diff_max, rho_ref_l2, tau_composite_hard, &
148 tau_composite_soft, tau_diff_l2, tau_diff_max, tau_ref_l2, volume
149 REAL(kind=dp), DIMENSION(3, 3) :: virial_xc_tmp
150 TYPE(cell_type), POINTER :: cell
151 TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: rho_ao_kp
152 TYPE(dft_control_type), POINTER :: dft_control
153 TYPE(mp_para_env_type), POINTER :: para_env
154 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
155 TYPE(pw_c1d_gs_type) :: tmp_g
156 TYPE(pw_c1d_gs_type), DIMENSION(:), POINTER :: rho_direct_g, rho_g, rho_hard_g, &
157 rho_m_gspace, rho_smooth_g, rho_smooth_model_g, rho_struct_g, tau_direct_g, tau_hard_g, &
158 tau_smooth_g, tau_smooth_model_g, tau_struct_g
159 TYPE(pw_c1d_gs_type), POINTER :: rho_nlcc_g, rho_nlcc_g_use, rho_nlcc_g_xc
160 TYPE(pw_env_type), POINTER :: pw_env
161 TYPE(pw_pool_type), POINTER :: auxbas_pw_pool, diagnostics_pw_pool, &
162 vdw_pw_pool, xc_pw_pool
163 TYPE(pw_r3d_rs_type), DIMENSION(:), POINTER :: my_vxc_rho, my_vxc_tau, rho_direct_r, &
164 rho_hard_r, rho_m_rspace, rho_r, rho_smooth_model_r, rho_smooth_r, rho_struct_r, tau, &
165 tau_direct_r, tau_hard_r, tau_smooth_model_r, tau_smooth_r, tau_struct_r
166 TYPE(pw_r3d_rs_type), POINTER :: rho_nlcc, rho_nlcc_use, rho_nlcc_xc, &
167 tmp_pw, weights, weights_use, &
168 weights_xc
169 TYPE(section_vals_type), POINTER :: gauxc_section
170 TYPE(virial_type), POINTER :: virial
171
172 CALL timeset(routinen, handle)
173
174 cpassert(.NOT. ASSOCIATED(vxc_rho))
175 cpassert(.NOT. ASSOCIATED(vxc_tau))
176 NULLIFY (dft_control, pw_env, auxbas_pw_pool, diagnostics_pw_pool, xc_pw_pool, vdw_pw_pool, &
177 cell, my_vxc_rho, &
178 tmp_pw, my_vxc_tau, rho_g, rho_r, tau, rho_m_rspace, &
179 rho_m_gspace, rho_nlcc, rho_nlcc_g, rho_nlcc_g_use, rho_nlcc_g_xc, &
180 rho_nlcc_use, rho_nlcc_xc, rho_struct_r, rho_struct_g, tau_struct_g, tau_struct_r, &
181 weights_use, weights_xc, particle_set, rho_ao_kp, rho_direct_g, rho_direct_r, &
182 rho_hard_g, rho_hard_r, rho_smooth_g, rho_smooth_model_g, rho_smooth_model_r, &
183 rho_smooth_r, tau_direct_g, tau_direct_r, tau_hard_g, tau_hard_r, tau_smooth_g, &
184 tau_smooth_model_g, tau_smooth_model_r, tau_smooth_r, gauxc_section)
185
186 exc = 0.0_dp
187 my_just_energy = .false.
188 IF (PRESENT(just_energy)) my_just_energy = just_energy
189 my_adiabatic_rescale_factor = 1.0_dp
190 do_adiabatic_rescaling = .false.
191 IF (PRESENT(adiabatic_rescale_factor)) THEN
192 my_adiabatic_rescale_factor = adiabatic_rescale_factor
193 do_adiabatic_rescaling = .true.
194 END IF
195
196 CALL get_ks_env(ks_env, &
197 dft_control=dft_control, &
198 pw_env=pw_env, &
199 cell=cell, &
200 para_env=para_env, &
201 particle_set=particle_set, &
202 xcint_weights=weights, &
203 virial=virial, &
204 rho_nlcc=rho_nlcc, &
205 rho_nlcc_g=rho_nlcc_g)
206 rho_nlcc_use => rho_nlcc
207 rho_nlcc_g_use => rho_nlcc_g
208 weights_use => weights
209
210 CALL qs_rho_get(rho_struct, &
211 tau_r_valid=tau_r_valid, &
212 tau_g_valid=tau_g_valid, &
213 rho_g_valid=rho_g_valid, &
214 rho_r=rho_struct_r, &
215 rho_g=rho_struct_g, &
216 tau_g=tau_struct_g, &
217 tau_r=tau_struct_r)
218 rho_smooth_r => rho_struct_r
219 rho_smooth_g => rho_struct_g
220 tau_smooth_r => tau_struct_r
221 tau_smooth_g => tau_struct_g
222
223 compute_virial = virial%pv_calculate .AND. (.NOT. virial%pv_numer)
224 IF (compute_virial) THEN
225 virial%pv_xc = 0.0_dp
226 END IF
227
228 CALL section_vals_val_get(xc_section, "XC_FUNCTIONAL%_SECTION_PARAMETERS_", &
229 i_val=myfun)
230 CALL section_vals_val_get(xc_section, "VDW_POTENTIAL%POTENTIAL_TYPE", &
231 i_val=vdw)
232
233 vdw_nl = (vdw == xc_vdw_fun_nonloc)
234 ! this combination has not been investigated
235 cpassert(.NOT. (do_adiabatic_rescaling .AND. vdw_nl))
236 ! are the necessary inputs available
237 IF (.NOT. (PRESENT(dispersion_env) .AND. PRESENT(edisp))) THEN
238 vdw_nl = .false.
239 END IF
240 IF (PRESENT(edisp)) edisp = 0.0_dp
241 native_skala_grid = xc_section_uses_native_skala_evaluator(xc_section)
242 defer_native_skala_to_atom_composite = .false.
243 IF (PRESENT(native_skala_defer_to_atom_composite)) THEN
244 defer_native_skala_to_atom_composite = native_skala_defer_to_atom_composite
245 END IF
246 cpassert(.NOT. defer_native_skala_to_atom_composite .OR. native_skala_grid)
247 native_gapw_composite_reference = native_skala_grid .AND. &
248 native_skala_gapw_composite_reference(xc_section) .AND. &
249 (dft_control%qs_control%gapw .OR. &
250 dft_control%qs_control%gapw_xc)
251 IF (PRESENT(native_gapw_composite_override)) THEN
252 native_gapw_composite_reference = native_skala_grid .AND. &
253 native_gapw_composite_override .AND. &
254 (dft_control%qs_control%gapw .OR. &
255 dft_control%qs_control%gapw_xc)
256 END IF
257 native_gapw_composite_direct_ao = native_gapw_composite_reference .AND. &
258 native_skala_gapw_composite_direct_ao(xc_section)
259 native_grid_diagnostics = .false.
260 IF (native_gapw_composite_reference) THEN
261 gauxc_section => get_gauxc_section(xc_section)
262 cpassert(ASSOCIATED(gauxc_section))
263 CALL section_vals_val_get(gauxc_section, "NATIVE_GRID_DIAGNOSTICS", &
264 l_val=native_grid_diagnostics)
265 END IF
266 IF (native_gapw_composite_reference) NULLIFY (weights_use)
267
268 IF (myfun /= xc_none .OR. vdw_nl) THEN
269
270 ! test if the real space density is available
271 cpassert(ASSOCIATED(rho_struct))
272 IF (dft_control%nspins /= 1 .AND. dft_control%nspins /= 2) THEN
273 cpabort("nspins must be 1 or 2")
274 END IF
275 mspin = SIZE(rho_struct_r)
276 IF (dft_control%nspins == 2 .AND. mspin == 1) THEN
277 cpabort("Spin count mismatch")
278 END IF
279
280 ! there are some options related to SIC here.
281 ! Normal DFT computes E(rho_alpha,rho_beta) (or its variant E(2*rho_alpha) for non-LSD)
282 ! SIC can E(rho_alpha,rho_beta)-b*(E(rho_alpha,rho_beta)-E(rho_beta,rho_beta))
283 ! or compute E(rho_alpha,rho_beta)-b*E(rho_alpha-rho_beta,0)
284
285 ! my_scaling is the scaling needed of the standard E(rho_alpha,rho_beta) term
286 my_scaling = 1.0_dp
287 SELECT CASE (dft_control%sic_method_id)
288 CASE (sic_none)
289 ! all fine
290 CASE (sic_mauri_spz, sic_ad)
291 ! no idea yet what to do here in that case
292 cpassert(.NOT. tau_r_valid)
293 CASE (sic_mauri_us)
294 my_scaling = 1.0_dp - dft_control%sic_scaling_b
295 ! no idea yet what to do here in that case
296 cpassert(.NOT. tau_r_valid)
297 CASE (sic_eo)
298 ! NOTHING TO BE DONE
299 CASE DEFAULT
300 ! this case has not yet been treated here
301 cpabort("NYI")
302 END SELECT
303
304 IF (dft_control%sic_scaling_b == 0.0_dp) THEN
305 sic_scaling_b_zero = .true.
306 ELSE
307 sic_scaling_b_zero = .false.
308 END IF
309
310 IF (PRESENT(pw_env_external)) THEN
311 pw_env => pw_env_external
312 END IF
313 CALL pw_env_get(pw_env, xc_pw_pool=xc_pw_pool, auxbas_pw_pool=auxbas_pw_pool)
314
315 IF (native_gapw_composite_reference) THEN
316 cpassert(tau_r_valid)
317 cpassert(PRESENT(qs_env_external))
318 ALLOCATE (rho_hard_r(mspin), rho_hard_g(mspin), tau_hard_r(mspin), tau_hard_g(mspin))
319 DO ispin = 1, mspin
320 CALL auxbas_pw_pool%create_pw(rho_hard_r(ispin))
321 CALL auxbas_pw_pool%create_pw(rho_hard_g(ispin))
322 CALL auxbas_pw_pool%create_pw(tau_hard_r(ispin))
323 CALL auxbas_pw_pool%create_pw(tau_hard_g(ispin))
324 END DO
325 IF (native_gapw_composite_direct_ao) THEN
326 CALL qs_rho_get(rho_struct, rho_ao_kp=rho_ao_kp)
327 cpassert(ASSOCIATED(rho_ao_kp))
328 cpassert(SIZE(rho_ao_kp, 2) == 1)
329 DO ispin = 1, mspin
330 CALL calculate_rho_elec( &
331 matrix_p=rho_ao_kp(ispin, 1)%matrix, rho=rho_hard_r(ispin), &
332 rho_gspace=rho_hard_g(ispin), ks_env=ks_env)
333 CALL calculate_rho_elec( &
334 matrix_p=rho_ao_kp(ispin, 1)%matrix, rho=tau_hard_r(ispin), &
335 rho_gspace=tau_hard_g(ispin), ks_env=ks_env, compute_tau=.true.)
336 END DO
337 ELSE
338 CALL get_pw_grid_info(rho_hard_g(1)%pw_grid, vol=volume)
339 q_max = sqrt(maxval(rho_hard_g(1)%pw_grid%gsq))
340 CALL calculate_rhotot_elec_gspace( &
341 qs_env=qs_env_external, auxbas_pw_pool=auxbas_pw_pool, &
342 rhotot_elec_gspace=rho_hard_g(1), q_max=q_max, &
343 rho_hard=composite_hard_integral, rho_soft=composite_soft_integral, &
344 rho_source=rho_struct, allow_nonorthorhombic=.true.)
345 rho_composite_hard = composite_hard_integral
346 rho_composite_soft = composite_soft_integral
347 CALL calculate_rhotot_elec_gspace( &
348 qs_env=qs_env_external, auxbas_pw_pool=auxbas_pw_pool, &
349 rhotot_elec_gspace=tau_hard_g(1), q_max=q_max, &
350 rho_hard=composite_hard_integral, rho_soft=composite_soft_integral, &
351 compute_tau=.true., rho_source=rho_struct, allow_nonorthorhombic=.true.)
352 tau_composite_hard = composite_hard_integral
353 tau_composite_soft = composite_soft_integral
354 IF (para_env%mepos == 0) THEN
355 output_unit = cp_logger_get_default_io_unit()
356 IF (output_unit > 0) THEN
357 WRITE (unit=output_unit, fmt="(/,T2,A,2(1X,ES19.11))") &
358 "SKALA_GPW| Composite rho hard and soft integrals", &
359 rho_composite_hard, rho_composite_soft
360 WRITE (unit=output_unit, fmt="(T2,A,2(1X,ES19.11))") &
361 "SKALA_GPW| Composite tau hard and soft integrals", &
362 tau_composite_hard, tau_composite_soft
363 END IF
364 END IF
365 IF (mspin == 1) THEN
366 CALL pw_scale(rho_hard_g(1), 1.0_dp/volume)
367 CALL pw_scale(tau_hard_g(1), 1.0_dp/volume)
368 ELSE
369 CALL calculate_rhotot_elec_gspace( &
370 qs_env=qs_env_external, auxbas_pw_pool=auxbas_pw_pool, &
371 rhotot_elec_gspace=rho_hard_g(2), q_max=q_max, &
372 rho_hard=composite_hard_integral, rho_soft=composite_soft_integral, &
373 fsign=-1.0_dp, rho_source=rho_struct, allow_nonorthorhombic=.true.)
374 CALL calculate_rhotot_elec_gspace( &
375 qs_env=qs_env_external, auxbas_pw_pool=auxbas_pw_pool, &
376 rhotot_elec_gspace=tau_hard_g(2), q_max=q_max, &
377 rho_hard=composite_hard_integral, rho_soft=composite_soft_integral, &
378 fsign=-1.0_dp, compute_tau=.true., rho_source=rho_struct, &
379 allow_nonorthorhombic=.true.)
380 CALL pw_scale(rho_hard_g(1), 0.5_dp/volume)
381 CALL pw_scale(rho_hard_g(2), 0.5_dp/volume)
382 CALL auxbas_pw_pool%create_pw(tmp_g)
383 CALL pw_copy(rho_hard_g(1), tmp_g)
384 CALL pw_axpy(rho_hard_g(2), rho_hard_g(1), 1.0_dp)
385 CALL pw_axpy(rho_hard_g(2), tmp_g, -1.0_dp)
386 CALL pw_copy(tmp_g, rho_hard_g(2))
387 CALL pw_scale(tau_hard_g(1), 0.5_dp/volume)
388 CALL pw_scale(tau_hard_g(2), 0.5_dp/volume)
389 CALL pw_copy(tau_hard_g(1), tmp_g)
390 CALL pw_axpy(tau_hard_g(2), tau_hard_g(1), 1.0_dp)
391 CALL pw_axpy(tau_hard_g(2), tmp_g, -1.0_dp)
392 CALL pw_copy(tmp_g, tau_hard_g(2))
393 CALL auxbas_pw_pool%give_back_pw(tmp_g)
394 END IF
395 END IF
396 DO ispin = 1, mspin
397 CALL pw_transfer(rho_hard_g(ispin), rho_hard_r(ispin))
398 CALL pw_transfer(tau_hard_g(ispin), tau_hard_r(ispin))
399 END DO
400 IF (native_grid_diagnostics .AND. .NOT. native_gapw_composite_direct_ao) THEN
401 CALL qs_rho_get(rho_struct, rho_ao_kp=rho_ao_kp)
402 cpassert(ASSOCIATED(rho_ao_kp))
403 cpassert(SIZE(rho_ao_kp, 2) == 1)
404 ALLOCATE (rho_direct_r(mspin), rho_direct_g(mspin), &
405 tau_direct_r(mspin), tau_direct_g(mspin))
406 rho_diff_l2 = 0.0_dp
407 rho_diff_max = 0.0_dp
408 rho_ref_l2 = 0.0_dp
409 tau_diff_l2 = 0.0_dp
410 tau_diff_max = 0.0_dp
411 tau_ref_l2 = 0.0_dp
412 direct_rho_integral = 0.0_dp
413 direct_tau_integral = 0.0_dp
414 DO ispin = 1, mspin
415 CALL auxbas_pw_pool%create_pw(rho_direct_r(ispin))
416 CALL auxbas_pw_pool%create_pw(rho_direct_g(ispin))
417 CALL auxbas_pw_pool%create_pw(tau_direct_r(ispin))
418 CALL auxbas_pw_pool%create_pw(tau_direct_g(ispin))
419 CALL calculate_rho_elec( &
420 matrix_p=rho_ao_kp(ispin, 1)%matrix, rho=rho_direct_r(ispin), &
421 rho_gspace=rho_direct_g(ispin), ks_env=ks_env)
422 CALL calculate_rho_elec( &
423 matrix_p=rho_ao_kp(ispin, 1)%matrix, rho=tau_direct_r(ispin), &
424 rho_gspace=tau_direct_g(ispin), ks_env=ks_env, compute_tau=.true.)
425 direct_rho_integral = direct_rho_integral + &
426 pw_integrate_function(rho_direct_r(ispin))
427 direct_tau_integral = direct_tau_integral + &
428 pw_integrate_function(tau_direct_r(ispin))
429 DO k = lbound(rho_direct_r(ispin)%array, 3), ubound(rho_direct_r(ispin)%array, 3)
430 DO j = lbound(rho_direct_r(ispin)%array, 2), ubound(rho_direct_r(ispin)%array, 2)
431 DO i = lbound(rho_direct_r(ispin)%array, 1), ubound(rho_direct_r(ispin)%array, 1)
432 delta = rho_hard_r(ispin)%array(i, j, k) - &
433 rho_direct_r(ispin)%array(i, j, k)
434 rho_diff_l2 = rho_diff_l2 + delta*delta
435 rho_diff_max = max(rho_diff_max, abs(delta))
436 rho_ref_l2 = rho_ref_l2 + rho_direct_r(ispin)%array(i, j, k)**2
437 delta = tau_hard_r(ispin)%array(i, j, k) - &
438 tau_direct_r(ispin)%array(i, j, k)
439 tau_diff_l2 = tau_diff_l2 + delta*delta
440 tau_diff_max = max(tau_diff_max, abs(delta))
441 tau_ref_l2 = tau_ref_l2 + tau_direct_r(ispin)%array(i, j, k)**2
442 END DO
443 END DO
444 END DO
445 END DO
446 CALL para_env%sum(rho_diff_l2)
447 CALL para_env%sum(rho_ref_l2)
448 CALL para_env%sum(tau_diff_l2)
449 CALL para_env%sum(tau_ref_l2)
450 CALL para_env%max(rho_diff_max)
451 CALL para_env%max(tau_diff_max)
452 rho_diff_l2 = sqrt(rho_diff_l2/max(rho_ref_l2, tiny(1.0_dp)))
453 tau_diff_l2 = sqrt(tau_diff_l2/max(tau_ref_l2, tiny(1.0_dp)))
454 IF (para_env%mepos == 0) THEN
455 output_unit = cp_logger_get_default_io_unit()
456 IF (output_unit > 0) THEN
457 WRITE (unit=output_unit, fmt="(T2,A,2(1X,ES19.11))") &
458 "SKALA_GPW| Direct AO rho and tau integrals", &
459 direct_rho_integral, direct_tau_integral
460 WRITE (unit=output_unit, fmt="(T2,A,2(1X,ES19.11))") &
461 "SKALA_GPW| Composite/direct relative L2 rho and tau", &
462 rho_diff_l2, tau_diff_l2
463 WRITE (unit=output_unit, fmt="(T2,A,2(1X,ES19.11))") &
464 "SKALA_GPW| Composite/direct Linf rho and tau", &
465 rho_diff_max, tau_diff_max
466 END IF
467 END IF
468 DO ispin = 1, mspin
469 CALL auxbas_pw_pool%give_back_pw(rho_direct_r(ispin))
470 CALL auxbas_pw_pool%give_back_pw(rho_direct_g(ispin))
471 CALL auxbas_pw_pool%give_back_pw(tau_direct_r(ispin))
472 CALL auxbas_pw_pool%give_back_pw(tau_direct_g(ispin))
473 END DO
474 DEALLOCATE (rho_direct_r, rho_direct_g, tau_direct_r, tau_direct_g)
475 END IF
476 composite_rho_r_integral = 0.0_dp
477 composite_tau_r_integral = 0.0_dp
478 composite_rho_min = huge(1.0_dp)
479 composite_rho_max = -huge(1.0_dp)
480 composite_tau_min = huge(1.0_dp)
481 composite_tau_max = -huge(1.0_dp)
482 DO ispin = 1, mspin
483 composite_rho_r_integral = composite_rho_r_integral + &
484 pw_integrate_function(rho_hard_r(ispin))
485 composite_tau_r_integral = composite_tau_r_integral + &
486 pw_integrate_function(tau_hard_r(ispin))
487 composite_rho_min = min(composite_rho_min, minval(rho_hard_r(ispin)%array))
488 composite_rho_max = max(composite_rho_max, maxval(rho_hard_r(ispin)%array))
489 composite_tau_min = min(composite_tau_min, minval(tau_hard_r(ispin)%array))
490 composite_tau_max = max(composite_tau_max, maxval(tau_hard_r(ispin)%array))
491 END DO
492 CALL para_env%min(composite_rho_min)
493 CALL para_env%max(composite_rho_max)
494 CALL para_env%min(composite_tau_min)
495 CALL para_env%max(composite_tau_max)
496 IF (para_env%mepos == 0) THEN
497 output_unit = cp_logger_get_default_io_unit()
498 IF (output_unit > 0) THEN
499 WRITE (unit=output_unit, fmt="(T2,A,1X,ES19.11)") &
500 "SKALA_GPW| Composite real-grid rho integral", composite_rho_r_integral
501 WRITE (unit=output_unit, fmt="(T2,A,1X,ES19.11)") &
502 "SKALA_GPW| Composite real-grid tau integral", composite_tau_r_integral
503 WRITE (unit=output_unit, fmt="(T2,A,2(1X,ES19.11))") &
504 "SKALA_GPW| Composite real-grid rho min and max", &
505 composite_rho_min, composite_rho_max
506 WRITE (unit=output_unit, fmt="(T2,A,2(1X,ES19.11))") &
507 "SKALA_GPW| Composite real-grid tau min and max", &
508 composite_tau_min, composite_tau_max
509 END IF
510 END IF
511 rho_struct_r => rho_hard_r
512 rho_struct_g => rho_hard_g
513 tau_struct_r => tau_hard_r
514 tau_struct_g => tau_hard_g
515 rho_g_valid = .true.
516 tau_g_valid = .true.
517 END IF
518
519 uf_grid = .NOT. pw_grid_compare(auxbas_pw_pool%pw_grid, xc_pw_pool%pw_grid)
520
521 IF (.NOT. uf_grid) THEN
522 rho_r => rho_struct_r
523
524 IF (tau_r_valid) THEN
525 tau => tau_struct_r
526 END IF
527
528 ! for gradient corrected functional the density in g space might
529 ! be useful so if we have it, we pass it in
530 IF (rho_g_valid) THEN
531 rho_g => rho_struct_g
532 END IF
533 ELSE
534 cpassert(rho_g_valid)
535 ALLOCATE (rho_r(mspin))
536 ALLOCATE (rho_g(mspin))
537 DO ispin = 1, mspin
538 CALL xc_pw_pool%create_pw(rho_g(ispin))
539 CALL pw_transfer(rho_struct_g(ispin), rho_g(ispin))
540 END DO
541 DO ispin = 1, mspin
542 CALL xc_pw_pool%create_pw(rho_r(ispin))
543 CALL pw_transfer(rho_g(ispin), rho_r(ispin))
544 END DO
545 IF (tau_r_valid) THEN
546 ALLOCATE (tau(mspin))
547 DO ispin = 1, mspin
548 CALL xc_pw_pool%create_pw(tau(ispin))
549 block
550 TYPE(pw_c1d_gs_type) :: tau_g_aux, tau_g_xc
551 CALL xc_pw_pool%create_pw(tau_g_xc)
552 IF (tau_g_valid) THEN
553 CALL pw_transfer(tau_struct_g(ispin), tau_g_xc)
554 ELSE
555 CALL auxbas_pw_pool%create_pw(tau_g_aux)
556 CALL pw_transfer(tau_struct_r(ispin), tau_g_aux)
557 CALL pw_transfer(tau_g_aux, tau_g_xc)
558 CALL auxbas_pw_pool%give_back_pw(tau_g_aux)
559 END IF
560 CALL pw_transfer(tau_g_xc, tau(ispin))
561 CALL xc_pw_pool%give_back_pw(tau_g_xc)
562 END block
563 END DO
564 END IF
565 IF (ASSOCIATED(weights) .AND. .NOT. native_gapw_composite_reference) THEN
566 ALLOCATE (weights_xc)
567 CALL xc_pw_pool%create_pw(weights_xc)
568 block
569 TYPE(pw_c1d_gs_type) :: weights_g_aux, weights_g_xc
570 CALL auxbas_pw_pool%create_pw(weights_g_aux)
571 CALL xc_pw_pool%create_pw(weights_g_xc)
572 CALL pw_transfer(weights, weights_g_aux)
573 CALL pw_transfer(weights_g_aux, weights_g_xc)
574 CALL pw_transfer(weights_g_xc, weights_xc)
575 CALL xc_pw_pool%give_back_pw(weights_g_xc)
576 CALL auxbas_pw_pool%give_back_pw(weights_g_aux)
577 END block
578 weights_use => weights_xc
579 END IF
580 IF (ASSOCIATED(rho_nlcc)) THEN
581 cpassert(ASSOCIATED(rho_nlcc_g))
582 ALLOCATE (rho_nlcc_g_xc, rho_nlcc_xc)
583 CALL xc_pw_pool%create_pw(rho_nlcc_g_xc)
584 CALL xc_pw_pool%create_pw(rho_nlcc_xc)
585 CALL pw_transfer(rho_nlcc_g, rho_nlcc_g_xc)
586 CALL pw_transfer(rho_nlcc_g_xc, rho_nlcc_xc)
587 rho_nlcc_use => rho_nlcc_xc
588 rho_nlcc_g_use => rho_nlcc_g_xc
589 END IF
590 END IF
591
592 IF (native_gapw_composite_reference) THEN
593 composite_rho_r_integral = 0.0_dp
594 composite_tau_r_integral = 0.0_dp
595 DO ispin = 1, mspin
596 composite_rho_r_integral = composite_rho_r_integral + &
597 pw_integrate_function(rho_r(ispin))
598 composite_tau_r_integral = composite_tau_r_integral + &
599 pw_integrate_function(tau(ispin))
600 END DO
601 IF (para_env%mepos == 0) THEN
602 output_unit = cp_logger_get_default_io_unit()
603 IF (output_unit > 0) THEN
604 WRITE (unit=output_unit, fmt="(T2,A,1X,ES19.11)") &
605 "SKALA_GPW| Composite XC-grid rho integral", composite_rho_r_integral
606 WRITE (unit=output_unit, fmt="(T2,A,1X,ES19.11)") &
607 "SKALA_GPW| Composite XC-grid tau integral", composite_tau_r_integral
608 END IF
609 END IF
610 END IF
611
612 ! NLCC augments rho and grad(rho); tau remains valence-only.
613 IF (ASSOCIATED(rho_nlcc_use)) THEN
614 factor = 1.0_dp
615 DO ispin = 1, mspin
616 CALL pw_axpy(rho_nlcc_use, rho_r(ispin), factor)
617 CALL pw_axpy(rho_nlcc_g_use, rho_g(ispin), factor)
618 END DO
619 END IF
620
621 !
622 ! here the rho_r, rho_g, tau is what it should be
623 ! we get back the right my_vxc_rho and my_vxc_tau as required
624 !
625 IF (defer_native_skala_to_atom_composite) THEN
626 ! The atom-centered composite path supplies the complete SKALA energy and
627 ! VXC below qs_vxc_create. Keep only storage for that VXC (and any additive
628 ! nonlocal-dispersion potential) instead of evaluating the smooth grid here.
629 ALLOCATE (my_vxc_rho(mspin), my_vxc_tau(mspin))
630 DO ispin = 1, mspin
631 CALL xc_pw_pool%create_pw(my_vxc_rho(ispin))
632 CALL xc_pw_pool%create_pw(my_vxc_tau(ispin))
633 CALL pw_zero(my_vxc_rho(ispin))
634 CALL pw_zero(my_vxc_tau(ispin))
635 END DO
636 ELSE IF (native_skala_grid) THEN
637 CALL skala_gpw_eval(vxc_rho=my_vxc_rho, vxc_tau=my_vxc_tau, exc=exc, &
638 rho_r=rho_r, rho_g=rho_g, tau=tau, xc_section=xc_section, &
639 weights=weights_use, pw_pool=xc_pw_pool, &
640 particle_set=particle_set, cell=cell, &
641 compute_virial=compute_virial, virial_xc=virial%pv_xc, &
642 just_energy=my_just_energy, atom_force=native_skala_atom_force)
643 ELSE IF (my_just_energy) THEN
644 exc = xc_exc_calc(rho_r=rho_r, tau=tau, &
645 rho_g=rho_g, xc_section=xc_section, &
646 weights=weights_use, pw_pool=xc_pw_pool)
647
648 ELSE
649 CALL xc_vxc_pw_create(vxc_rho=my_vxc_rho, vxc_tau=my_vxc_tau, rho_r=rho_r, &
650 rho_g=rho_g, tau=tau, exc=exc, &
651 xc_section=xc_section, &
652 weights=weights_use, pw_pool=xc_pw_pool, &
653 compute_virial=compute_virial, &
654 virial_xc=virial%pv_xc)
655 END IF
656
657 ! The composite density is assembled through an R-to-G-to-R projection. Apply the
658 ! discrete adjoint of that projection before integrating VXC with the smooth AO density.
659 IF (native_gapw_composite_reference .AND. .NOT. my_just_energy) THEN
660 IF (uf_grid) THEN
661 diagnostics_pw_pool => xc_pw_pool
662 ELSE
663 diagnostics_pw_pool => auxbas_pw_pool
664 END IF
665 CALL diagnostics_pw_pool%create_pw(tmp_g)
666 DO ispin = 1, mspin
667 CALL pw_transfer(my_vxc_rho(ispin), tmp_g)
668 CALL pw_transfer(tmp_g, my_vxc_rho(ispin))
669 CALL pw_transfer(my_vxc_tau(ispin), tmp_g)
670 CALL pw_transfer(tmp_g, my_vxc_tau(ispin))
671 END DO
672 CALL diagnostics_pw_pool%give_back_pw(tmp_g)
673 NULLIFY (diagnostics_pw_pool)
674 END IF
675
676 IF (native_gapw_composite_reference .AND. native_grid_diagnostics .AND. &
677 .NOT. my_just_energy) THEN
678 IF (uf_grid) THEN
679 diagnostics_pw_pool => xc_pw_pool
680 ELSE
681 diagnostics_pw_pool => auxbas_pw_pool
682 END IF
683 ALLOCATE (rho_smooth_model_g(mspin), rho_smooth_model_r(mspin), &
684 tau_smooth_model_g(mspin), tau_smooth_model_r(mspin))
685 DO ispin = 1, mspin
686 CALL diagnostics_pw_pool%create_pw(rho_smooth_model_g(ispin))
687 CALL diagnostics_pw_pool%create_pw(rho_smooth_model_r(ispin))
688 CALL diagnostics_pw_pool%create_pw(tau_smooth_model_g(ispin))
689 CALL diagnostics_pw_pool%create_pw(tau_smooth_model_r(ispin))
690 IF (uf_grid) THEN
691 IF (rho_g_valid) THEN
692 CALL pw_transfer(rho_smooth_g(ispin), rho_smooth_model_g(ispin))
693 ELSE
694 CALL auxbas_pw_pool%create_pw(tmp_g)
695 CALL pw_transfer(rho_smooth_r(ispin), tmp_g)
696 CALL pw_transfer(tmp_g, rho_smooth_model_g(ispin))
697 CALL auxbas_pw_pool%give_back_pw(tmp_g)
698 END IF
699 IF (tau_g_valid) THEN
700 CALL pw_transfer(tau_smooth_g(ispin), tau_smooth_model_g(ispin))
701 ELSE
702 CALL auxbas_pw_pool%create_pw(tmp_g)
703 CALL pw_transfer(tau_smooth_r(ispin), tmp_g)
704 CALL pw_transfer(tmp_g, tau_smooth_model_g(ispin))
705 CALL auxbas_pw_pool%give_back_pw(tmp_g)
706 END IF
707 ELSE
708 IF (rho_g_valid) THEN
709 CALL pw_copy(rho_smooth_g(ispin), rho_smooth_model_g(ispin))
710 ELSE
711 CALL pw_transfer(rho_smooth_r(ispin), rho_smooth_model_g(ispin))
712 END IF
713 IF (tau_g_valid) THEN
714 CALL pw_copy(tau_smooth_g(ispin), tau_smooth_model_g(ispin))
715 ELSE
716 CALL pw_transfer(tau_smooth_r(ispin), tau_smooth_model_g(ispin))
717 END IF
718 END IF
719 CALL pw_transfer(rho_smooth_model_g(ispin), rho_smooth_model_r(ispin))
720 CALL pw_transfer(tau_smooth_model_g(ispin), tau_smooth_model_r(ispin))
721 IF (ASSOCIATED(rho_nlcc_use)) THEN
722 CALL pw_axpy(rho_nlcc_use, rho_smooth_model_r(ispin), 1.0_dp)
723 CALL pw_axpy(rho_nlcc_g_use, rho_smooth_model_g(ispin), 1.0_dp)
724 END IF
725 END DO
726 CALL diagnose_gapw_composite_direction( &
727 rho_r, rho_g, tau, rho_smooth_model_r, rho_smooth_model_g, &
728 tau_smooth_model_r, my_vxc_rho, my_vxc_tau, xc_section, weights_use, &
729 diagnostics_pw_pool, particle_set, cell, para_env)
730 DO ispin = 1, mspin
731 CALL diagnostics_pw_pool%give_back_pw(rho_smooth_model_g(ispin))
732 CALL diagnostics_pw_pool%give_back_pw(rho_smooth_model_r(ispin))
733 CALL diagnostics_pw_pool%give_back_pw(tau_smooth_model_g(ispin))
734 CALL diagnostics_pw_pool%give_back_pw(tau_smooth_model_r(ispin))
735 END DO
736 DEALLOCATE (rho_smooth_model_g, rho_smooth_model_r, &
737 tau_smooth_model_g, tau_smooth_model_r)
738 NULLIFY (diagnostics_pw_pool)
739 END IF
740
741 ! remove the nlcc densities (keep stuff in original state)
742 IF (ASSOCIATED(rho_nlcc_use)) THEN
743 factor = -1.0_dp
744 DO ispin = 1, mspin
745 CALL pw_axpy(rho_nlcc_use, rho_r(ispin), factor)
746 CALL pw_axpy(rho_nlcc_g_use, rho_g(ispin), factor)
747 END DO
748 END IF
749
750 ! calclulate non-local vdW functional
751 ! only if this XC_SECTION has it
752 ! if yes, we use the dispersion_env from ks_env
753 ! this is dangerous, as it assumes a special connection xc_section -> qs_env
754 IF (vdw_nl) THEN
755 CALL get_ks_env(ks_env=ks_env, para_env=para_env)
756 ! no SIC functionals allowed
757 cpassert(dft_control%sic_method_id == sic_none)
758 !
759 CALL pw_env_get(pw_env, vdw_pw_pool=vdw_pw_pool)
760 IF (my_just_energy) THEN
761 CALL calculate_dispersion_nonloc(my_vxc_rho, rho_r, rho_g, edisp, dispersion_env, &
762 my_just_energy, vdw_pw_pool, xc_pw_pool, para_env)
763 ELSE
764 CALL calculate_dispersion_nonloc(my_vxc_rho, rho_r, rho_g, edisp, dispersion_env, &
765 my_just_energy, vdw_pw_pool, xc_pw_pool, para_env, virial=virial)
766 END IF
767 END IF
768
769 !! Apply rescaling to the potential if requested
770 IF (.NOT. my_just_energy) THEN
771 IF (do_adiabatic_rescaling) THEN
772 IF (ASSOCIATED(my_vxc_rho)) THEN
773 DO ispin = 1, SIZE(my_vxc_rho)
774 CALL pw_scale(my_vxc_rho(ispin), my_adiabatic_rescale_factor)
775 END DO
776 END IF
777 END IF
778 END IF
779
780 IF (my_scaling /= 1.0_dp) THEN
781 exc = exc*my_scaling
782 IF (ASSOCIATED(my_vxc_rho)) THEN
783 DO ispin = 1, SIZE(my_vxc_rho)
784 CALL pw_scale(my_vxc_rho(ispin), my_scaling)
785 END DO
786 END IF
787 IF (ASSOCIATED(my_vxc_tau)) THEN
788 DO ispin = 1, SIZE(my_vxc_tau)
789 CALL pw_scale(my_vxc_tau(ispin), my_scaling)
790 END DO
791 END IF
792 END IF
793
794 ! we have pw data for the xc, qs_ks requests coeff structure, here we transfer
795 ! pw -> coeff
796 IF (ASSOCIATED(my_vxc_rho)) THEN
797 vxc_rho => my_vxc_rho
798 NULLIFY (my_vxc_rho)
799 END IF
800 IF (ASSOCIATED(my_vxc_tau)) THEN
801 vxc_tau => my_vxc_tau
802 NULLIFY (my_vxc_tau)
803 END IF
804 IF (uf_grid) THEN
805 DO ispin = 1, SIZE(rho_r)
806 CALL xc_pw_pool%give_back_pw(rho_r(ispin))
807 END DO
808 DEALLOCATE (rho_r)
809 IF (ASSOCIATED(rho_g)) THEN
810 DO ispin = 1, SIZE(rho_g)
811 CALL xc_pw_pool%give_back_pw(rho_g(ispin))
812 END DO
813 DEALLOCATE (rho_g)
814 END IF
815 END IF
816
817 ! compute again the xc but now for Exc(m,o) and the opposite sign
818 IF (dft_control%sic_method_id == sic_mauri_spz .AND. .NOT. sic_scaling_b_zero) THEN
819 ALLOCATE (rho_m_rspace(2), rho_m_gspace(2))
820 CALL xc_pw_pool%create_pw(rho_m_gspace(1))
821 CALL xc_pw_pool%create_pw(rho_m_rspace(1))
822 CALL pw_copy(rho_struct_r(1), rho_m_rspace(1))
823 CALL pw_axpy(rho_struct_r(2), rho_m_rspace(1), alpha=-1._dp)
824 CALL pw_copy(rho_struct_g(1), rho_m_gspace(1))
825 CALL pw_axpy(rho_struct_g(2), rho_m_gspace(1), alpha=-1._dp)
826 ! bit sad, these will be just zero...
827 CALL xc_pw_pool%create_pw(rho_m_gspace(2))
828 CALL xc_pw_pool%create_pw(rho_m_rspace(2))
829 CALL pw_zero(rho_m_rspace(2))
830 CALL pw_zero(rho_m_gspace(2))
831
832 IF (my_just_energy) THEN
833 exc_m = xc_exc_calc(rho_r=rho_m_rspace, tau=tau, &
834 rho_g=rho_m_gspace, xc_section=xc_section, &
835 weights=weights_use, pw_pool=xc_pw_pool)
836 ELSE
837 ! virial untested
838 cpassert(.NOT. compute_virial)
839 CALL xc_vxc_pw_create(vxc_rho=my_vxc_rho, vxc_tau=my_vxc_tau, rho_r=rho_m_rspace, &
840 rho_g=rho_m_gspace, tau=tau, exc=exc_m, &
841 xc_section=xc_section, &
842 weights=weights_use, pw_pool=xc_pw_pool, &
843 compute_virial=.false., &
844 virial_xc=virial_xc_tmp)
845 END IF
846
847 exc = exc - dft_control%sic_scaling_b*exc_m
848
849 ! and take care of the potential only vxc_rho is taken into account
850 IF (.NOT. my_just_energy) THEN
851 CALL pw_axpy(my_vxc_rho(1), vxc_rho(1), -dft_control%sic_scaling_b)
852 CALL pw_axpy(my_vxc_rho(1), vxc_rho(2), dft_control%sic_scaling_b)
853 CALL my_vxc_rho(1)%release()
854 CALL my_vxc_rho(2)%release()
855 DEALLOCATE (my_vxc_rho)
856 END IF
857
858 DO ispin = 1, 2
859 CALL xc_pw_pool%give_back_pw(rho_m_rspace(ispin))
860 CALL xc_pw_pool%give_back_pw(rho_m_gspace(ispin))
861 END DO
862 DEALLOCATE (rho_m_rspace)
863 DEALLOCATE (rho_m_gspace)
864
865 END IF
866
867 ! now we have - sum_s N_s * Exc(rho_s/N_s,0)
868 IF (dft_control%sic_method_id == sic_ad .AND. .NOT. sic_scaling_b_zero) THEN
869
870 ! find out how many elecs we have
871 CALL get_ks_env(ks_env, nelectron_spin=nelec_spin)
872
873 ALLOCATE (rho_m_rspace(2), rho_m_gspace(2))
874 DO ispin = 1, 2
875 CALL xc_pw_pool%create_pw(rho_m_gspace(ispin))
876 CALL xc_pw_pool%create_pw(rho_m_rspace(ispin))
877 END DO
878
879 DO ispin = 1, 2
880 IF (nelec_spin(ispin) > 0.0_dp) THEN
881 nelec_s_inv = 1.0_dp/nelec_spin(ispin)
882 ELSE
883 ! does it matter if there are no electrons with this spin (H) ?
884 nelec_s_inv = 0.0_dp
885 END IF
886 CALL pw_copy(rho_struct_r(ispin), rho_m_rspace(1))
887 CALL pw_copy(rho_struct_g(ispin), rho_m_gspace(1))
888 CALL pw_scale(rho_m_rspace(1), nelec_s_inv)
889 CALL pw_scale(rho_m_gspace(1), nelec_s_inv)
890 CALL pw_zero(rho_m_rspace(2))
891 CALL pw_zero(rho_m_gspace(2))
892
893 IF (my_just_energy) THEN
894 exc_m = xc_exc_calc(rho_r=rho_m_rspace, tau=tau, &
895 rho_g=rho_m_gspace, xc_section=xc_section, &
896 weights=weights_use, pw_pool=xc_pw_pool)
897 ELSE
898 ! virial untested
899 cpassert(.NOT. compute_virial)
900 CALL xc_vxc_pw_create(vxc_rho=my_vxc_rho, vxc_tau=my_vxc_tau, rho_r=rho_m_rspace, &
901 rho_g=rho_m_gspace, tau=tau, exc=exc_m, &
902 xc_section=xc_section, &
903 weights=weights_use, pw_pool=xc_pw_pool, &
904 compute_virial=.false., &
905 virial_xc=virial_xc_tmp)
906 END IF
907
908 exc = exc - dft_control%sic_scaling_b*nelec_spin(ispin)*exc_m
909
910 ! and take care of the potential only vxc_rho is taken into account
911 IF (.NOT. my_just_energy) THEN
912 CALL pw_axpy(my_vxc_rho(1), vxc_rho(ispin), -dft_control%sic_scaling_b)
913 CALL my_vxc_rho(1)%release()
914 CALL my_vxc_rho(2)%release()
915 DEALLOCATE (my_vxc_rho)
916 END IF
917 END DO
918
919 DO ispin = 1, 2
920 CALL xc_pw_pool%give_back_pw(rho_m_rspace(ispin))
921 CALL xc_pw_pool%give_back_pw(rho_m_gspace(ispin))
922 END DO
923 DEALLOCATE (rho_m_rspace)
924 DEALLOCATE (rho_m_gspace)
925
926 END IF
927
928 ! compute again the xc but now for Exc(n_down,n_down)
929 IF (dft_control%sic_method_id == sic_mauri_us .AND. .NOT. sic_scaling_b_zero) THEN
930 ALLOCATE (rho_r(2))
931 rho_r(1) = rho_struct_r(2)
932 rho_r(2) = rho_struct_r(2)
933 IF (rho_g_valid) THEN
934 ALLOCATE (rho_g(2))
935 rho_g(1) = rho_struct_g(2)
936 rho_g(2) = rho_struct_g(2)
937 END IF
938
939 IF (my_just_energy) THEN
940 exc_m = xc_exc_calc(rho_r=rho_r, tau=tau, &
941 rho_g=rho_g, xc_section=xc_section, &
942 weights=weights_use, pw_pool=xc_pw_pool)
943 ELSE
944 ! virial untested
945 cpassert(.NOT. compute_virial)
946 CALL xc_vxc_pw_create(vxc_rho=my_vxc_rho, vxc_tau=my_vxc_tau, rho_r=rho_r, &
947 rho_g=rho_g, tau=tau, exc=exc_m, &
948 xc_section=xc_section, &
949 weights=weights_use, pw_pool=xc_pw_pool, &
950 compute_virial=.false., &
951 virial_xc=virial_xc_tmp)
952 END IF
953
954 exc = exc + dft_control%sic_scaling_b*exc_m
955
956 ! and take care of the potential
957 IF (.NOT. my_just_energy) THEN
958 ! both go to minority spin
959 CALL pw_axpy(my_vxc_rho(1), vxc_rho(2), 2.0_dp*dft_control%sic_scaling_b)
960 CALL my_vxc_rho(1)%release()
961 CALL my_vxc_rho(2)%release()
962 DEALLOCATE (my_vxc_rho)
963 END IF
964 DEALLOCATE (rho_r, rho_g)
965
966 END IF
967
968 !
969 ! cleanups
970 !
971 IF (uf_grid .AND. (ASSOCIATED(vxc_rho) .OR. ASSOCIATED(vxc_tau))) THEN
972 block
973 TYPE(pw_r3d_rs_type) :: tmp_pw
974 TYPE(pw_c1d_gs_type) :: tmp_g, tmp_g2
975 CALL xc_pw_pool%create_pw(tmp_g)
976 CALL auxbas_pw_pool%create_pw(tmp_g2)
977 IF (ASSOCIATED(vxc_rho)) THEN
978 DO ispin = 1, SIZE(vxc_rho)
979 CALL auxbas_pw_pool%create_pw(tmp_pw)
980 CALL pw_transfer(vxc_rho(ispin), tmp_g)
981 CALL pw_transfer(tmp_g, tmp_g2)
982 CALL pw_transfer(tmp_g2, tmp_pw)
983 CALL xc_pw_pool%give_back_pw(vxc_rho(ispin))
984 vxc_rho(ispin) = tmp_pw
985 END DO
986 END IF
987 IF (ASSOCIATED(vxc_tau)) THEN
988 DO ispin = 1, SIZE(vxc_tau)
989 CALL auxbas_pw_pool%create_pw(tmp_pw)
990 CALL pw_transfer(vxc_tau(ispin), tmp_g)
991 CALL pw_transfer(tmp_g, tmp_g2)
992 CALL pw_transfer(tmp_g2, tmp_pw)
993 CALL xc_pw_pool%give_back_pw(vxc_tau(ispin))
994 vxc_tau(ispin) = tmp_pw
995 END DO
996 END IF
997 CALL auxbas_pw_pool%give_back_pw(tmp_g2)
998 CALL xc_pw_pool%give_back_pw(tmp_g)
999 END block
1000 END IF
1001 IF (ASSOCIATED(tau) .AND. uf_grid) THEN
1002 DO ispin = 1, SIZE(tau)
1003 CALL xc_pw_pool%give_back_pw(tau(ispin))
1004 END DO
1005 DEALLOCATE (tau)
1006 END IF
1007 IF (ASSOCIATED(weights_xc)) THEN
1008 CALL xc_pw_pool%give_back_pw(weights_xc)
1009 DEALLOCATE (weights_xc)
1010 END IF
1011 IF (ASSOCIATED(rho_nlcc_xc)) THEN
1012 CALL xc_pw_pool%give_back_pw(rho_nlcc_xc)
1013 DEALLOCATE (rho_nlcc_xc)
1014 END IF
1015 IF (ASSOCIATED(rho_nlcc_g_xc)) THEN
1016 CALL xc_pw_pool%give_back_pw(rho_nlcc_g_xc)
1017 DEALLOCATE (rho_nlcc_g_xc)
1018 END IF
1019
1020 IF (ASSOCIATED(rho_hard_r)) THEN
1021 DO ispin = 1, SIZE(rho_hard_r)
1022 CALL auxbas_pw_pool%give_back_pw(rho_hard_r(ispin))
1023 CALL auxbas_pw_pool%give_back_pw(rho_hard_g(ispin))
1024 CALL auxbas_pw_pool%give_back_pw(tau_hard_r(ispin))
1025 CALL auxbas_pw_pool%give_back_pw(tau_hard_g(ispin))
1026 END DO
1027 DEALLOCATE (rho_hard_r, rho_hard_g, tau_hard_r, tau_hard_g)
1028 END IF
1029
1030 END IF
1031
1032 CALL timestop(handle)
1033
1034 END SUBROUTINE qs_vxc_create
1035
1036! **************************************************************************************************
1037!> \brief Check the SKALA energy derivative along the reconstructed GAPW one-center correction.
1038!> \param rho_combined_r combined smooth plus hard-minus-soft density
1039!> \param rho_combined_g reciprocal representation of rho_combined_r
1040!> \param tau_combined combined kinetic-energy density
1041!> \param rho_smooth_r smooth density on the model grid
1042!> \param rho_smooth_g reciprocal representation of rho_smooth_r
1043!> \param tau_smooth smooth kinetic-energy density on the model grid
1044!> \param vxc_rho derivative of the SKALA energy with respect to rho
1045!> \param vxc_tau derivative of the SKALA energy with respect to tau
1046!> \param xc_section ...
1047!> \param weights ...
1048!> \param pw_pool ...
1049!> \param particle_set ...
1050!> \param cell ...
1051!> \param para_env ...
1052! **************************************************************************************************
1053 SUBROUTINE diagnose_gapw_composite_direction( &
1054 rho_combined_r, rho_combined_g, tau_combined, rho_smooth_r, rho_smooth_g, tau_smooth, &
1055 vxc_rho, vxc_tau, xc_section, weights, pw_pool, particle_set, cell, para_env)
1056
1057 TYPE(pw_r3d_rs_type), DIMENSION(:), POINTER :: rho_combined_r
1058 TYPE(pw_c1d_gs_type), DIMENSION(:), POINTER :: rho_combined_g
1059 TYPE(pw_r3d_rs_type), DIMENSION(:), POINTER :: tau_combined, rho_smooth_r
1060 TYPE(pw_c1d_gs_type), DIMENSION(:), POINTER :: rho_smooth_g
1061 TYPE(pw_r3d_rs_type), DIMENSION(:), POINTER :: tau_smooth, vxc_rho, vxc_tau
1062 TYPE(section_vals_type), POINTER :: xc_section
1063 TYPE(pw_r3d_rs_type), POINTER :: weights
1064 TYPE(pw_pool_type), POINTER :: pw_pool
1065 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
1066 TYPE(cell_type), POINTER :: cell
1067 TYPE(mp_para_env_type), POINTER :: para_env
1068
1069 REAL(kind=dp), DIMENSION(3), PARAMETER :: fd_steps = [1.0e-2_dp, 3.0e-3_dp, 1.0e-3_dp]
1070
1071 INTEGER :: ispin, istep, nspin, output_unit
1072 REAL(kind=dp) :: exc_minus, exc_plus, fd_derivative, fd_step, rho_contraction, &
1073 rho_g_contraction, tau_contraction, tau_g_contraction, total_contraction
1074 REAL(kind=dp), DIMENSION(3, 3) :: virial_dummy
1075 TYPE(pw_c1d_gs_type), DIMENSION(:), POINTER :: rho_test_g, tau_test_g, vxc_rho_g, &
1076 vxc_tau_g
1077 TYPE(pw_r3d_rs_type), DIMENSION(:), POINTER :: rho_test_r, tau_test, vxc_dummy_rho, &
1078 vxc_dummy_tau
1079
1080 NULLIFY (rho_test_g, rho_test_r, tau_test, tau_test_g, vxc_dummy_rho, vxc_dummy_tau, &
1081 vxc_rho_g, vxc_tau_g)
1082 nspin = SIZE(rho_combined_r)
1083 cpassert(SIZE(rho_combined_g) == nspin)
1084 cpassert(SIZE(tau_combined) == nspin)
1085 cpassert(SIZE(rho_smooth_r) == nspin)
1086 cpassert(SIZE(rho_smooth_g) == nspin)
1087 cpassert(SIZE(tau_smooth) == nspin)
1088 cpassert(SIZE(vxc_rho) == nspin)
1089 cpassert(SIZE(vxc_tau) == nspin)
1090
1091 ALLOCATE (rho_test_g(nspin), rho_test_r(nspin), tau_test(nspin), tau_test_g(nspin), &
1092 vxc_rho_g(nspin), vxc_tau_g(nspin))
1093 DO ispin = 1, nspin
1094 CALL pw_pool%create_pw(rho_test_g(ispin))
1095 CALL pw_pool%create_pw(rho_test_r(ispin))
1096 CALL pw_pool%create_pw(tau_test(ispin))
1097 CALL pw_pool%create_pw(tau_test_g(ispin))
1098 CALL pw_pool%create_pw(vxc_rho_g(ispin))
1099 CALL pw_pool%create_pw(vxc_tau_g(ispin))
1100 CALL pw_copy(rho_combined_g(ispin), rho_test_g(ispin))
1101 CALL pw_axpy(rho_smooth_g(ispin), rho_test_g(ispin), -1.0_dp)
1102 CALL pw_copy(rho_combined_r(ispin), rho_test_r(ispin))
1103 CALL pw_axpy(rho_smooth_r(ispin), rho_test_r(ispin), -1.0_dp)
1104 CALL pw_copy(tau_combined(ispin), tau_test(ispin))
1105 CALL pw_axpy(tau_smooth(ispin), tau_test(ispin), -1.0_dp)
1106 CALL pw_transfer(tau_test(ispin), tau_test_g(ispin))
1107 CALL pw_transfer(vxc_rho(ispin), vxc_rho_g(ispin))
1108 CALL pw_transfer(vxc_tau(ispin), vxc_tau_g(ispin))
1109 END DO
1110 rho_contraction = 0.0_dp
1111 rho_g_contraction = 0.0_dp
1112 tau_contraction = 0.0_dp
1113 tau_g_contraction = 0.0_dp
1114 DO ispin = 1, nspin
1115 rho_contraction = rho_contraction + pw_integral_ab(rho_test_r(ispin), vxc_rho(ispin))
1116 rho_g_contraction = rho_g_contraction + pw_integral_ab(rho_test_g(ispin), vxc_rho_g(ispin))
1117 tau_contraction = tau_contraction + pw_integral_ab(tau_test(ispin), vxc_tau(ispin))
1118 tau_g_contraction = tau_g_contraction + pw_integral_ab(tau_test_g(ispin), vxc_tau_g(ispin))
1119 END DO
1120
1121 total_contraction = rho_contraction + tau_contraction
1122 DO istep = 1, SIZE(fd_steps)
1123 fd_step = fd_steps(istep)
1124 DO ispin = 1, nspin
1125 CALL pw_copy(rho_combined_g(ispin), rho_test_g(ispin))
1126 CALL pw_axpy(rho_combined_g(ispin), rho_test_g(ispin), fd_step)
1127 CALL pw_axpy(rho_smooth_g(ispin), rho_test_g(ispin), -fd_step)
1128 CALL pw_copy(rho_combined_r(ispin), rho_test_r(ispin))
1129 CALL pw_axpy(rho_combined_r(ispin), rho_test_r(ispin), fd_step)
1130 CALL pw_axpy(rho_smooth_r(ispin), rho_test_r(ispin), -fd_step)
1131 CALL pw_copy(tau_combined(ispin), tau_test(ispin))
1132 CALL pw_axpy(tau_combined(ispin), tau_test(ispin), fd_step)
1133 CALL pw_axpy(tau_smooth(ispin), tau_test(ispin), -fd_step)
1134 END DO
1135 CALL skala_gpw_eval( &
1136 vxc_rho=vxc_dummy_rho, vxc_tau=vxc_dummy_tau, exc=exc_plus, rho_r=rho_test_r, &
1137 rho_g=rho_test_g, tau=tau_test, xc_section=xc_section, weights=weights, pw_pool=pw_pool, &
1138 particle_set=particle_set, cell=cell, compute_virial=.false., virial_xc=virial_dummy, &
1139 just_energy=.true.)
1140 cpassert(.NOT. ASSOCIATED(vxc_dummy_rho))
1141 cpassert(.NOT. ASSOCIATED(vxc_dummy_tau))
1142
1143 DO ispin = 1, nspin
1144 CALL pw_copy(rho_combined_g(ispin), rho_test_g(ispin))
1145 CALL pw_axpy(rho_combined_g(ispin), rho_test_g(ispin), -fd_step)
1146 CALL pw_axpy(rho_smooth_g(ispin), rho_test_g(ispin), fd_step)
1147 CALL pw_copy(rho_combined_r(ispin), rho_test_r(ispin))
1148 CALL pw_axpy(rho_combined_r(ispin), rho_test_r(ispin), -fd_step)
1149 CALL pw_axpy(rho_smooth_r(ispin), rho_test_r(ispin), fd_step)
1150 CALL pw_copy(tau_combined(ispin), tau_test(ispin))
1151 CALL pw_axpy(tau_combined(ispin), tau_test(ispin), -fd_step)
1152 CALL pw_axpy(tau_smooth(ispin), tau_test(ispin), fd_step)
1153 END DO
1154 CALL skala_gpw_eval( &
1155 vxc_rho=vxc_dummy_rho, vxc_tau=vxc_dummy_tau, exc=exc_minus, rho_r=rho_test_r, &
1156 rho_g=rho_test_g, tau=tau_test, xc_section=xc_section, weights=weights, pw_pool=pw_pool, &
1157 particle_set=particle_set, cell=cell, compute_virial=.false., virial_xc=virial_dummy, &
1158 just_energy=.true.)
1159 cpassert(.NOT. ASSOCIATED(vxc_dummy_rho))
1160 cpassert(.NOT. ASSOCIATED(vxc_dummy_tau))
1161
1162 fd_derivative = (exc_plus - exc_minus)/(2.0_dp*fd_step)
1163 IF (para_env%mepos == 0) THEN
1164 output_unit = cp_logger_get_default_io_unit()
1165 IF (output_unit > 0) THEN
1166 WRITE (output_unit, "(/,T2,A,1X,ES12.4,2(1X,ES20.10))") &
1167 "SKALA_GPW| Composite-direction FD step and energies", &
1168 fd_step, exc_plus, exc_minus
1169 WRITE (output_unit, "(T2,A,1X,ES20.10)") &
1170 "SKALA_GPW| Composite-direction energy FD", fd_derivative
1171 WRITE (output_unit, "(T2,A,1X,ES20.10)") &
1172 "SKALA_GPW| Composite-direction FD/VXC difference", &
1173 fd_derivative - total_contraction
1174 END IF
1175 END IF
1176 END DO
1177 IF (para_env%mepos == 0) THEN
1178 output_unit = cp_logger_get_default_io_unit()
1179 IF (output_unit > 0) THEN
1180 WRITE (output_unit, "(T2,A,1X,ES20.10)") &
1181 "SKALA_GPW| Composite-direction rho VXC contraction", rho_contraction
1182 WRITE (output_unit, "(T2,A,1X,ES20.10)") &
1183 "SKALA_GPW| Composite-direction rho G-space contraction", rho_g_contraction
1184 WRITE (output_unit, "(T2,A,1X,ES20.10)") &
1185 "SKALA_GPW| Composite-direction tau VXC contraction", tau_contraction
1186 WRITE (output_unit, "(T2,A,1X,ES20.10)") &
1187 "SKALA_GPW| Composite-direction tau G-space contraction", tau_g_contraction
1188 WRITE (output_unit, "(T2,A,1X,ES20.10)") &
1189 "SKALA_GPW| Composite-direction total VXC contraction", total_contraction
1190 END IF
1191 END IF
1192
1193 DO ispin = 1, nspin
1194 CALL pw_pool%give_back_pw(rho_test_g(ispin))
1195 CALL pw_pool%give_back_pw(rho_test_r(ispin))
1196 CALL pw_pool%give_back_pw(tau_test(ispin))
1197 CALL pw_pool%give_back_pw(tau_test_g(ispin))
1198 CALL pw_pool%give_back_pw(vxc_rho_g(ispin))
1199 CALL pw_pool%give_back_pw(vxc_tau_g(ispin))
1200 END DO
1201 DEALLOCATE (rho_test_g, rho_test_r, tau_test, tau_test_g, vxc_rho_g, vxc_tau_g)
1202
1203 END SUBROUTINE diagnose_gapw_composite_direction
1204
1205! **************************************************************************************************
1206!> \brief calculates the XC density: E_xc(r) - V_xc(r)*rho(r) or E_xc(r)/rho(r)
1207!> \param ks_env to get all the needed things
1208!> \param rho_struct density
1209!> \param xc_section ...
1210!> \param dispersion_env ...
1211!> \param xc_ener will contain the xc energy density E_xc(r) - V_xc(r)*rho(r)
1212!> \param xc_den will contain the xc energy density E_xc(r)/rho(r)
1213!> \param exc will contain the xc energy density E_xc(r)
1214!> \param vxc ...
1215!> \param vtau ...
1216!> \author JGH
1217! **************************************************************************************************
1218 SUBROUTINE qs_xc_density(ks_env, rho_struct, xc_section, dispersion_env, &
1219 xc_ener, xc_den, exc, vxc, vtau)
1220
1221 TYPE(qs_ks_env_type), POINTER :: ks_env
1222 TYPE(qs_rho_type), POINTER :: rho_struct
1223 TYPE(section_vals_type), POINTER :: xc_section
1224 TYPE(qs_dispersion_type), OPTIONAL, POINTER :: dispersion_env
1225 TYPE(pw_r3d_rs_type), INTENT(INOUT), OPTIONAL :: xc_ener, xc_den
1226 TYPE(pw_r3d_rs_type), OPTIONAL :: exc
1227 TYPE(pw_r3d_rs_type), DIMENSION(:), OPTIONAL :: vxc, vtau
1228
1229 CHARACTER(len=*), PARAMETER :: routinen = 'qs_xc_density'
1230
1231 INTEGER :: handle, ispin, mspin, myfun, nspins, vdw
1232 LOGICAL :: rho_g_valid, tau_g_valid, tau_r_valid, &
1233 uf_grid, vdw_nl
1234 REAL(kind=dp) :: edisp, excint, factor, rho_cutoff
1235 REAL(kind=dp), DIMENSION(3, 3) :: vdum
1236 TYPE(cell_type), POINTER :: cell
1237 TYPE(dft_control_type), POINTER :: dft_control
1238 TYPE(mp_para_env_type), POINTER :: para_env
1239 TYPE(pw_c1d_gs_type), DIMENSION(:), POINTER :: rho_g, rho_struct_g, tau_g, tau_struct_g
1240 TYPE(pw_c1d_gs_type), POINTER :: rho_nlcc_g, rho_nlcc_g_use, rho_nlcc_g_xc
1241 TYPE(pw_env_type), POINTER :: pw_env
1242 TYPE(pw_pool_type), POINTER :: auxbas_pw_pool, vdw_pw_pool, xc_pw_pool
1243 TYPE(pw_r3d_rs_type) :: exc_r
1244 TYPE(pw_r3d_rs_type), DIMENSION(:), POINTER :: rho_r, rho_struct_r, tau_r, &
1245 tau_struct_r, vxc_rho, vxc_tau
1246 TYPE(pw_r3d_rs_type), POINTER :: rho_nlcc, rho_nlcc_use, rho_nlcc_xc, &
1247 weights, weights_use, weights_xc
1248
1249 CALL timeset(routinen, handle)
1250
1251 NULLIFY (dft_control, pw_env, auxbas_pw_pool, xc_pw_pool, vdw_pw_pool, cell, &
1252 rho_g, rho_struct_g, tau_g, tau_struct_g, rho_nlcc, rho_nlcc_g, &
1253 rho_nlcc_g_use, rho_nlcc_g_xc, rho_nlcc_use, rho_nlcc_xc, rho_r, &
1254 rho_struct_r, tau_r, tau_struct_r, vxc_rho, vxc_tau, weights, &
1255 weights_use, weights_xc)
1256
1257 CALL get_ks_env(ks_env, &
1258 dft_control=dft_control, &
1259 pw_env=pw_env, &
1260 cell=cell, &
1261 xcint_weights=weights, &
1262 rho_nlcc=rho_nlcc, &
1263 rho_nlcc_g=rho_nlcc_g)
1264
1265 CALL qs_rho_get(rho_struct, &
1266 tau_r_valid=tau_r_valid, &
1267 tau_g_valid=tau_g_valid, &
1268 rho_g_valid=rho_g_valid, &
1269 rho_r=rho_struct_r, &
1270 rho_g=rho_struct_g, &
1271 tau_r=tau_struct_r, &
1272 tau_g=tau_struct_g)
1273 nspins = dft_control%nspins
1274 mspin = SIZE(rho_struct_r)
1275 rho_r => rho_struct_r
1276 rho_g => rho_struct_g
1277 tau_r => tau_struct_r
1278 tau_g => tau_struct_g
1279 rho_nlcc_use => rho_nlcc
1280 rho_nlcc_g_use => rho_nlcc_g
1281 weights_use => weights
1282
1283 CALL section_vals_val_get(xc_section, "XC_FUNCTIONAL%_SECTION_PARAMETERS_", i_val=myfun)
1284 CALL section_vals_val_get(xc_section, "VDW_POTENTIAL%POTENTIAL_TYPE", i_val=vdw)
1285 vdw_nl = (vdw == xc_vdw_fun_nonloc)
1286 IF (PRESENT(xc_ener)) THEN
1287 IF (tau_r_valid) THEN
1288 CALL cp_warn(__location__, "Tau contribution will not be correctly handled")
1289 END IF
1290 END IF
1291 IF (vdw_nl) THEN
1292 CALL cp_warn(__location__, "vdW functional contribution will be ignored")
1293 END IF
1294
1295 CALL pw_env_get(pw_env, xc_pw_pool=xc_pw_pool, auxbas_pw_pool=auxbas_pw_pool)
1296 uf_grid = .NOT. pw_grid_compare(auxbas_pw_pool%pw_grid, xc_pw_pool%pw_grid)
1297
1298 IF (PRESENT(xc_ener)) THEN
1299 CALL pw_zero(xc_ener)
1300 END IF
1301 IF (PRESENT(xc_den)) THEN
1302 CALL pw_zero(xc_den)
1303 END IF
1304 IF (PRESENT(exc)) THEN
1305 CALL pw_zero(exc)
1306 END IF
1307 IF (PRESENT(vxc)) THEN
1308 DO ispin = 1, nspins
1309 CALL pw_zero(vxc(ispin))
1310 END DO
1311 END IF
1312 IF (PRESENT(vtau)) THEN
1313 DO ispin = 1, nspins
1314 CALL pw_zero(vtau(ispin))
1315 END DO
1316 END IF
1317
1318 IF (myfun /= xc_none) THEN
1319
1320 cpassert(ASSOCIATED(rho_struct))
1321 cpassert(dft_control%sic_method_id == sic_none)
1322
1323 IF (uf_grid) THEN
1324 NULLIFY (rho_r, rho_g, tau_r, tau_g)
1325 IF (rho_g_valid) THEN
1326 CALL create_density_on_pool(xc_pw_pool, rho_struct_g, rho_r, rho_g)
1327 ELSE IF (ASSOCIATED(rho_struct_r)) THEN
1328 CALL create_density_on_pool_from_r(auxbas_pw_pool, xc_pw_pool, rho_struct_r, rho_r, rho_g)
1329 ELSE
1330 cpabort("Fine Grid in qs_xc_density requires rho_r or rho_g")
1331 END IF
1332 IF (tau_r_valid) THEN
1333 IF (tau_g_valid) THEN
1334 CALL create_density_on_pool(xc_pw_pool, tau_struct_g, tau_r, tau_g)
1335 ELSE IF (ASSOCIATED(tau_struct_r)) THEN
1336 CALL create_density_on_pool_from_r(auxbas_pw_pool, xc_pw_pool, tau_struct_r, tau_r, tau_g)
1337 ELSE
1338 cpabort("Fine Grid in qs_xc_density requires tau_r or tau_g")
1339 END IF
1340 END IF
1341 IF (ASSOCIATED(weights)) THEN
1342 ALLOCATE (weights_xc)
1343 CALL xc_pw_pool%create_pw(weights_xc)
1344 CALL transfer_rspace_between_pools(auxbas_pw_pool, xc_pw_pool, weights, weights_xc)
1345 weights_use => weights_xc
1346 END IF
1347 IF (ASSOCIATED(rho_nlcc)) THEN
1348 cpassert(ASSOCIATED(rho_nlcc_g))
1349 ALLOCATE (rho_nlcc_g_xc, rho_nlcc_xc)
1350 CALL xc_pw_pool%create_pw(rho_nlcc_g_xc)
1351 CALL xc_pw_pool%create_pw(rho_nlcc_xc)
1352 CALL pw_transfer(rho_nlcc_g, rho_nlcc_g_xc)
1353 CALL pw_transfer(rho_nlcc_g_xc, rho_nlcc_xc)
1354 rho_nlcc_use => rho_nlcc_xc
1355 rho_nlcc_g_use => rho_nlcc_g_xc
1356 END IF
1357 END IF
1358
1359 ! add the nlcc densities
1360 IF (ASSOCIATED(rho_nlcc_use)) THEN
1361 factor = 1.0_dp
1362 DO ispin = 1, mspin
1363 CALL pw_axpy(rho_nlcc_use, rho_r(ispin), factor)
1364 CALL pw_axpy(rho_nlcc_g_use, rho_g(ispin), factor)
1365 END DO
1366 END IF
1367 NULLIFY (vxc_rho, vxc_tau)
1368 CALL xc_vxc_pw_create(vxc_rho=vxc_rho, vxc_tau=vxc_tau, rho_r=rho_r, &
1369 rho_g=rho_g, tau=tau_r, exc=excint, &
1370 xc_section=xc_section, &
1371 weights=weights_use, pw_pool=xc_pw_pool, &
1372 compute_virial=.false., &
1373 virial_xc=vdum, &
1374 exc_r=exc_r)
1375 ! calclulate non-local vdW functional
1376 ! only if this XC_SECTION has it
1377 ! if yes, we use the dispersion_env from ks_env
1378 ! this is dangerous, as it assumes a special connection xc_section -> qs_env
1379 IF (vdw_nl) THEN
1380 CALL get_ks_env(ks_env=ks_env, para_env=para_env)
1381 ! no SIC functionals allowed
1382 cpassert(dft_control%sic_method_id == sic_none)
1383 !
1384 CALL pw_env_get(pw_env, vdw_pw_pool=vdw_pw_pool)
1385 CALL calculate_dispersion_nonloc(vxc_rho, rho_r, rho_g, edisp, dispersion_env, &
1386 .false., vdw_pw_pool, xc_pw_pool, para_env)
1387 END IF
1388
1389 ! remove the nlcc densities (keep stuff in original state)
1390 IF (ASSOCIATED(rho_nlcc_use)) THEN
1391 factor = -1.0_dp
1392 DO ispin = 1, mspin
1393 CALL pw_axpy(rho_nlcc_use, rho_r(ispin), factor)
1394 CALL pw_axpy(rho_nlcc_g_use, rho_g(ispin), factor)
1395 END DO
1396 END IF
1397 !
1398 IF (PRESENT(xc_den)) THEN
1399 rho_cutoff = 1.e-14_dp
1400 IF (uf_grid) THEN
1401 block
1402 TYPE(pw_r3d_rs_type) :: tmp_pw
1403 CALL xc_pw_pool%create_pw(tmp_pw)
1404 CALL pw_copy(exc_r, tmp_pw)
1405 CALL calc_xc_density(tmp_pw, rho_r, rho_cutoff)
1406 CALL transfer_rspace_between_pools(xc_pw_pool, auxbas_pw_pool, tmp_pw, xc_den)
1407 CALL xc_pw_pool%give_back_pw(tmp_pw)
1408 END block
1409 ELSE
1410 CALL pw_copy(exc_r, xc_den)
1411 CALL calc_xc_density(xc_den, rho_r, rho_cutoff)
1412 END IF
1413 END IF
1414 IF (PRESENT(xc_ener)) THEN
1415 IF (uf_grid) THEN
1416 block
1417 TYPE(pw_r3d_rs_type) :: tmp_pw
1418 CALL xc_pw_pool%create_pw(tmp_pw)
1419 CALL pw_copy(exc_r, tmp_pw)
1420 DO ispin = 1, nspins
1421 CALL pw_multiply(tmp_pw, vxc_rho(ispin), rho_r(ispin), alpha=-1.0_dp)
1422 END DO
1423 CALL transfer_rspace_between_pools(xc_pw_pool, auxbas_pw_pool, tmp_pw, xc_ener)
1424 CALL xc_pw_pool%give_back_pw(tmp_pw)
1425 END block
1426 ELSE
1427 CALL pw_copy(exc_r, xc_ener)
1428 DO ispin = 1, nspins
1429 CALL pw_multiply(xc_ener, vxc_rho(ispin), rho_r(ispin), alpha=-1.0_dp)
1430 END DO
1431 END IF
1432 END IF
1433 IF (PRESENT(exc)) THEN
1434 IF (uf_grid) THEN
1435 CALL transfer_rspace_between_pools(xc_pw_pool, auxbas_pw_pool, exc_r, exc)
1436 ELSE
1437 CALL pw_copy(exc_r, exc)
1438 END IF
1439 END IF
1440 IF (PRESENT(vxc)) THEN
1441 DO ispin = 1, nspins
1442 IF (uf_grid) THEN
1443 CALL transfer_rspace_between_pools(xc_pw_pool, auxbas_pw_pool, vxc_rho(ispin), vxc(ispin))
1444 ELSE
1445 CALL pw_copy(vxc_rho(ispin), vxc(ispin))
1446 END IF
1447 END DO
1448 END IF
1449 IF (PRESENT(vtau) .AND. ASSOCIATED(vxc_tau)) THEN
1450 DO ispin = 1, nspins
1451 IF (uf_grid) THEN
1452 CALL transfer_rspace_between_pools(xc_pw_pool, auxbas_pw_pool, vxc_tau(ispin), vtau(ispin))
1453 ELSE
1454 CALL pw_copy(vxc_tau(ispin), vtau(ispin))
1455 END IF
1456 END DO
1457 END IF
1458 ! remove arrays
1459 IF (ASSOCIATED(vxc_rho)) THEN
1460 DO ispin = 1, nspins
1461 CALL vxc_rho(ispin)%release()
1462 END DO
1463 DEALLOCATE (vxc_rho)
1464 END IF
1465 IF (ASSOCIATED(vxc_tau)) THEN
1466 DO ispin = 1, nspins
1467 CALL vxc_tau(ispin)%release()
1468 END DO
1469 DEALLOCATE (vxc_tau)
1470 END IF
1471 CALL exc_r%release()
1472 IF (uf_grid) THEN
1473 CALL give_back_density_on_pool(xc_pw_pool, rho_r, rho_g)
1474 IF (ASSOCIATED(tau_r)) CALL give_back_density_on_pool(xc_pw_pool, tau_r, tau_g)
1475 IF (ASSOCIATED(weights_xc)) THEN
1476 CALL xc_pw_pool%give_back_pw(weights_xc)
1477 DEALLOCATE (weights_xc)
1478 END IF
1479 IF (ASSOCIATED(rho_nlcc_xc)) THEN
1480 CALL xc_pw_pool%give_back_pw(rho_nlcc_xc)
1481 DEALLOCATE (rho_nlcc_xc)
1482 END IF
1483 IF (ASSOCIATED(rho_nlcc_g_xc)) THEN
1484 CALL xc_pw_pool%give_back_pw(rho_nlcc_g_xc)
1485 DEALLOCATE (rho_nlcc_g_xc)
1486 END IF
1487 END IF
1488 !
1489 END IF
1490
1491 CALL timestop(handle)
1492
1493 END SUBROUTINE qs_xc_density
1494
1495! **************************************************************************************************
1496!> \brief transfers an r-space PW between two pools and writes into an existing target PW
1497!> \param source_pw_pool ...
1498!> \param target_pw_pool ...
1499!> \param source ...
1500!> \param TARGET ...
1501! **************************************************************************************************
1502 SUBROUTINE transfer_rspace_between_pools(source_pw_pool, target_pw_pool, source, TARGET)
1503 TYPE(pw_pool_type), POINTER :: source_pw_pool, target_pw_pool
1504 TYPE(pw_r3d_rs_type), INTENT(INOUT) :: source, target
1505
1506 TYPE(pw_c1d_gs_type) :: source_g, target_g
1507
1508 cpassert(ASSOCIATED(source_pw_pool))
1509 cpassert(ASSOCIATED(target_pw_pool))
1510
1511 IF (pw_grid_compare(source_pw_pool%pw_grid, target_pw_pool%pw_grid)) THEN
1512 CALL pw_copy(source, TARGET)
1513 ELSE
1514 CALL source_pw_pool%create_pw(source_g)
1515 CALL target_pw_pool%create_pw(target_g)
1516 CALL pw_transfer(source, source_g)
1517 CALL pw_transfer(source_g, target_g)
1518 CALL pw_transfer(target_g, TARGET)
1519 CALL target_pw_pool%give_back_pw(target_g)
1520 CALL source_pw_pool%give_back_pw(source_g)
1521 END IF
1522
1523 END SUBROUTINE transfer_rspace_between_pools
1524
1525! **************************************************************************************************
1526!> \brief transfers a g-space density to a given PW pool and creates its r-space representation
1527!> \param pw_pool ...
1528!> \param rho_g_in ...
1529!> \param rho_r_out ...
1530!> \param rho_g_out ...
1531! **************************************************************************************************
1532 SUBROUTINE create_density_on_pool(pw_pool, rho_g_in, rho_r_out, rho_g_out)
1533 TYPE(pw_pool_type), POINTER :: pw_pool
1534 TYPE(pw_c1d_gs_type), DIMENSION(:), POINTER :: rho_g_in
1535 TYPE(pw_r3d_rs_type), DIMENSION(:), POINTER :: rho_r_out
1536 TYPE(pw_c1d_gs_type), DIMENSION(:), POINTER :: rho_g_out
1537
1538 INTEGER :: ispin, nspins
1539
1540 cpassert(ASSOCIATED(pw_pool))
1541 cpassert(ASSOCIATED(rho_g_in))
1542
1543 nspins = SIZE(rho_g_in)
1544 ALLOCATE (rho_r_out(nspins), rho_g_out(nspins))
1545 DO ispin = 1, nspins
1546 CALL pw_pool%create_pw(rho_g_out(ispin))
1547 CALL pw_pool%create_pw(rho_r_out(ispin))
1548 CALL pw_transfer(rho_g_in(ispin), rho_g_out(ispin))
1549 CALL pw_transfer(rho_g_out(ispin), rho_r_out(ispin))
1550 END DO
1551
1552 END SUBROUTINE create_density_on_pool
1553
1554! **************************************************************************************************
1555!> \brief transfers an r-space density to a given PW pool and creates its g-space representation
1556!> \param source_pw_pool ...
1557!> \param target_pw_pool ...
1558!> \param rho_r_in ...
1559!> \param rho_r_out ...
1560!> \param rho_g_out ...
1561! **************************************************************************************************
1562 SUBROUTINE create_density_on_pool_from_r(source_pw_pool, target_pw_pool, rho_r_in, rho_r_out, rho_g_out)
1563 TYPE(pw_pool_type), POINTER :: source_pw_pool, target_pw_pool
1564 TYPE(pw_r3d_rs_type), DIMENSION(:), POINTER :: rho_r_in, rho_r_out
1565 TYPE(pw_c1d_gs_type), DIMENSION(:), POINTER :: rho_g_out
1566
1567 INTEGER :: ispin, nspins
1568 TYPE(pw_c1d_gs_type) :: rho_g_in
1569
1570 cpassert(ASSOCIATED(source_pw_pool))
1571 cpassert(ASSOCIATED(target_pw_pool))
1572 cpassert(ASSOCIATED(rho_r_in))
1573
1574 nspins = SIZE(rho_r_in)
1575 ALLOCATE (rho_r_out(nspins), rho_g_out(nspins))
1576 DO ispin = 1, nspins
1577 CALL source_pw_pool%create_pw(rho_g_in)
1578 CALL target_pw_pool%create_pw(rho_g_out(ispin))
1579 CALL target_pw_pool%create_pw(rho_r_out(ispin))
1580 CALL pw_transfer(rho_r_in(ispin), rho_g_in)
1581 CALL pw_transfer(rho_g_in, rho_g_out(ispin))
1582 CALL pw_transfer(rho_g_out(ispin), rho_r_out(ispin))
1583 CALL source_pw_pool%give_back_pw(rho_g_in)
1584 END DO
1585
1586 END SUBROUTINE create_density_on_pool_from_r
1587
1588! **************************************************************************************************
1589!> \brief returns temporary density arrays to the given PW pool
1590!> \param pw_pool ...
1591!> \param rho_r ...
1592!> \param rho_g ...
1593! **************************************************************************************************
1594 SUBROUTINE give_back_density_on_pool(pw_pool, rho_r, rho_g)
1595 TYPE(pw_pool_type), POINTER :: pw_pool
1596 TYPE(pw_r3d_rs_type), DIMENSION(:), POINTER :: rho_r
1597 TYPE(pw_c1d_gs_type), DIMENSION(:), POINTER :: rho_g
1598
1599 INTEGER :: ispin
1600
1601 cpassert(ASSOCIATED(pw_pool))
1602
1603 IF (ASSOCIATED(rho_r)) THEN
1604 DO ispin = 1, SIZE(rho_r)
1605 CALL pw_pool%give_back_pw(rho_r(ispin))
1606 END DO
1607 DEALLOCATE (rho_r)
1608 END IF
1609 IF (ASSOCIATED(rho_g)) THEN
1610 DO ispin = 1, SIZE(rho_g)
1611 CALL pw_pool%give_back_pw(rho_g(ispin))
1612 END DO
1613 DEALLOCATE (rho_g)
1614 END IF
1615
1616 END SUBROUTINE give_back_density_on_pool
1617
1618END MODULE qs_vxc
pw_methods::pw_axpy
Definition pw_methods.F:165
pw_methods::pw_copy
Definition pw_methods.F:146
pw_methods::pw_integral_ab
Definition pw_methods.F:222
pw_methods::pw_integrate_function
Definition pw_methods.F:116
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
cp_dbcsr_api
Definition cp_dbcsr_api.F:8
cp_log_handling
various routines to log and control the output. The idea is that decisions about where to log should ...
Definition cp_log_handling.F:41
cp_log_handling::cp_logger_get_default_io_unit
integer function, public cp_logger_get_default_io_unit(logger)
returns the unit nr for the ionode (-1 on all other processors) skips as well checks if the procs cal...
Definition cp_log_handling.F:523
gapw_gspace_reconstruction
GAPW reciprocal-space reconstruction and its discrete adjoint.
Definition gapw_gspace_reconstruction.F:12
gapw_gspace_reconstruction::calculate_rhotot_elec_gspace
subroutine, public calculate_rhotot_elec_gspace(qs_env, auxbas_pw_pool, rhotot_elec_gspace, q_max, rho_hard, rho_soft, fsign, compute_tau, rho_source, allow_nonorthorhombic)
The total electronic density in reciprocal space (g-space) is calculated.
Definition gapw_gspace_reconstruction.F:88
input_constants
collects all constants needed in input so that they can be used without circular dependencies
Definition input_constants.F:18
input_constants::sic_mauri_spz
integer, parameter, public sic_mauri_spz
Definition input_constants.F:602
input_constants::xc_vdw_fun_nonloc
integer, parameter, public xc_vdw_fun_nonloc
Definition input_constants.F:629
input_constants::sic_eo
integer, parameter, public sic_eo
Definition input_constants.F:602
input_constants::sic_mauri_us
integer, parameter, public sic_mauri_us
Definition input_constants.F:602
input_constants::sic_none
integer, parameter, public sic_none
Definition input_constants.F:602
input_constants::xc_none
integer, parameter, public xc_none
Definition input_constants.F:585
input_constants::sic_ad
integer, parameter, public sic_ad
Definition input_constants.F:602
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:1052
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:27
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_grids::get_pw_grid_info
subroutine, public get_pw_grid_info(pw_grid, id_nr, mode, vol, dvol, npts, ngpts, ngpts_cut, dr, cutoff, orthorhombic, gvectors, gsquare)
Access to information stored in the pw_grid_type.
Definition pw_grids.F:185
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_collocate_density
Calculate the plane wave density by collocating the primitive Gaussian functions (pgf).
Definition qs_collocate_density.F:38
qs_collocate_density::calculate_rho_elec
subroutine, public calculate_rho_elec(matrix_p, matrix_p_kp, rho, rho_gspace, total_rho, ks_env, soft_valid, compute_tau, compute_grad, basis_type, der_type, idir, task_list_external, pw_env_external)
computes the density corresponding to a given density matrix on the grid
Definition qs_collocate_density.F:1720
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:151
qs_dispersion_types
Definition of disperson types for DFT calculations.
Definition qs_dispersion_types.F:12
qs_environment_types
Definition qs_environment_types.F:14
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_vhxc, 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:343
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, qs_env_external, native_gapw_composite_override, native_skala_defer_to_atom_composite)
calculates and allocates the xc potential, already reducing it to the dependence on rho and the one o...
Definition qs_vxc.F:118
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:1220
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:404
skala_gpw_functional::native_skala_gapw_composite_direct_ao
logical function, public native_skala_gapw_composite_direct_ao(xc_section)
Return true if the GAPW composite reference uses direct full-ORB collocation.
Definition skala_gpw_functional.F:258
skala_gpw_functional::xc_section_uses_native_skala_evaluator
logical function, public xc_section_uses_native_skala_evaluator(xc_section)
Return true when SKALA must be evaluated by the CP2K-native grid machinery.
Definition skala_gpw_functional.F:164
skala_gpw_functional::get_gauxc_section
type(section_vals_type) function, pointer, public get_gauxc_section(xc_section)
Return the first GAUXC functional subsection, if present.
Definition skala_gpw_functional.F:3026
skala_gpw_functional::native_skala_gapw_composite_reference
logical function, public native_skala_gapw_composite_reference(xc_section)
Return true if native SKALA should use the full GAPW ORB density on one common grid.
Definition skala_gpw_functional.F:179
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:792
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:475
xc::calc_xc_density
subroutine, public calc_xc_density(pot, rho, rho_cutoff)
Definition xc.F:403
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:758
cp_dbcsr_api::dbcsr_p_type
Definition cp_dbcsr_api.F:172
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:767
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_environment_types::qs_environment_type
Definition qs_environment_types.F:222
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