Line data Source code
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 module that builds the second order perturbation kernel
10 : !> kpp1 = delta_rho|_P delta_rho|_P E drho(P1) drho
11 : !> \par History
12 : !> 07.2002 created [fawzi]
13 : !> \author Fawzi Mohamed
14 : ! **************************************************************************************************
15 : MODULE qs_kpp1_env_methods
16 : USE admm_types, ONLY: admm_type,&
17 : get_admm_env
18 : USE atomic_kind_types, ONLY: atomic_kind_type
19 : USE cp_control_types, ONLY: dft_control_type
20 : USE cp_dbcsr_api, ONLY: dbcsr_add,&
21 : dbcsr_copy,&
22 : dbcsr_p_type,&
23 : dbcsr_set
24 : USE cp_dbcsr_operations, ONLY: dbcsr_allocate_matrix_set
25 : USE cp_log_handling, ONLY: cp_to_string
26 : USE hartree_local_methods, ONLY: Vh_1c_gg_integrals
27 : USE input_constants, ONLY: do_admm_aux_exch_func_none
28 : USE input_section_types, ONLY: section_vals_type
29 : USE kahan_sum, ONLY: accurate_sum
30 : USE kinds, ONLY: dp
31 : USE lri_environment_types, ONLY: lri_density_type,&
32 : lri_environment_type,&
33 : lri_kind_type
34 : USE lri_ks_methods, ONLY: calculate_lri_ks_matrix
35 : USE message_passing, ONLY: mp_para_env_type
36 : USE pw_env_types, ONLY: pw_env_get,&
37 : pw_env_type
38 : USE pw_methods, ONLY: pw_axpy,&
39 : pw_copy,&
40 : pw_integrate_function,&
41 : pw_scale,&
42 : pw_transfer,&
43 : pw_zero
44 : USE pw_poisson_methods, ONLY: pw_poisson_solve
45 : USE pw_poisson_types, ONLY: pw_poisson_type
46 : USE pw_pool_types, ONLY: pw_pool_type
47 : USE pw_types, ONLY: pw_c1d_gs_type,&
48 : pw_r3d_rs_type
49 : USE qs_dispersion_nonloc, ONLY: qs_fxc_nlvdw_apply
50 : USE qs_environment_types, ONLY: get_qs_env,&
51 : qs_environment_type
52 : USE qs_fxc, ONLY: qs_fxc_apply,&
53 : qs_fxc_prep
54 : USE qs_gapw_densities, ONLY: prepare_gapw_den
55 : USE qs_integrate_potential, ONLY: integrate_v_rspace,&
56 : integrate_v_rspace_diagonal,&
57 : integrate_v_rspace_one_center
58 : USE qs_kpp1_env_types, ONLY: qs_kpp1_env_type
59 : USE qs_ks_atom, ONLY: update_ks_atom
60 : USE qs_p_env_types, ONLY: qs_p_env_type
61 : USE qs_rho0_ggrid, ONLY: integrate_vhg0_rspace
62 : USE qs_rho_atom_types, ONLY: rho_atom_type
63 : USE qs_rho_types, ONLY: qs_rho_get,&
64 : qs_rho_type
65 : #include "./base/base_uses.f90"
66 :
67 : IMPLICIT NONE
68 :
69 : PRIVATE
70 :
71 : LOGICAL, PRIVATE, PARAMETER :: debug_this_module = .TRUE.
72 : CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'qs_kpp1_env_methods'
73 :
74 : PUBLIC :: kpp1_create, kpp1_check_i_alloc, calc_kpp1
75 :
76 : CONTAINS
77 :
78 : ! **************************************************************************************************
79 : !> \brief allocates and initializes a kpp1_env
80 : !> \param kpp1_env the environment to initialize
81 : !> \par History
82 : !> 07.2002 created [fawzi]
83 : !> \author Fawzi Mohamed
84 : ! **************************************************************************************************
85 1878 : SUBROUTINE kpp1_create(kpp1_env)
86 : TYPE(qs_kpp1_env_type) :: kpp1_env
87 :
88 1878 : NULLIFY (kpp1_env%v_ao)
89 :
90 1878 : END SUBROUTINE kpp1_create
91 :
92 : ! **************************************************************************************************
93 : !> \brief ...
94 : !> \param rho1_xc ...
95 : !> \param rho1 ...
96 : !> \param xc_section ...
97 : !> \param lrigpw ...
98 : !> \param qs_env ...
99 : !> \param p_env ...
100 : !> \param calc_forces ...
101 : !> \param calc_virial ...
102 : !> \param virial ...
103 : ! **************************************************************************************************
104 1794 : SUBROUTINE calc_kpp1(rho1_xc, rho1, xc_section, lrigpw, qs_env, p_env, &
105 : calc_forces, calc_virial, virial)
106 :
107 : TYPE(qs_rho_type), POINTER :: rho1_xc, rho1
108 : TYPE(section_vals_type), POINTER :: xc_section
109 : LOGICAL, INTENT(IN) :: lrigpw
110 : TYPE(qs_environment_type), POINTER :: qs_env
111 : TYPE(qs_p_env_type) :: p_env
112 : LOGICAL, INTENT(IN), OPTIONAL :: calc_forces, calc_virial
113 : REAL(KIND=dp), DIMENSION(3, 3), INTENT(INOUT), &
114 : OPTIONAL :: virial
115 :
116 : CHARACTER(len=*), PARAMETER :: routineN = 'calc_kpp1'
117 :
118 : INTEGER :: handle, ikind, ispin, nkind, ns, nspins
119 : LOGICAL :: do_onecenter, gapw, gapw_xc, &
120 : my_calc_forces
121 : REAL(KIND=dp) :: alpha, energy_hartree, energy_hartree_1c
122 1794 : TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
123 1794 : TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: k1mat, rho_ao
124 1794 : TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: ksmat, psmat
125 : TYPE(dft_control_type), POINTER :: dft_control
126 : TYPE(lri_density_type), POINTER :: lri_density
127 : TYPE(lri_environment_type), POINTER :: lri_env
128 1794 : TYPE(lri_kind_type), DIMENSION(:), POINTER :: lri_v_int
129 : TYPE(mp_para_env_type), POINTER :: para_env
130 : TYPE(pw_c1d_gs_type) :: rho1_tot_gspace
131 1794 : TYPE(pw_c1d_gs_type), DIMENSION(:), POINTER :: rho1_g
132 : TYPE(pw_env_type), POINTER :: pw_env
133 : TYPE(pw_poisson_type), POINTER :: poisson_env
134 : TYPE(pw_pool_type), POINTER :: auxbas_pw_pool
135 : TYPE(pw_r3d_rs_type) :: v_hartree_rspace
136 1794 : TYPE(pw_r3d_rs_type), DIMENSION(:), POINTER :: rho0_r, rho1_r, v_rspace_new, v_xc, &
137 1794 : v_xc_tau
138 : TYPE(qs_kpp1_env_type), POINTER :: kpp1_env
139 : TYPE(qs_rho_type), POINTER :: rho, rho0_fxc, rho1_fxc
140 1794 : TYPE(rho_atom_type), DIMENSION(:), POINTER :: rho1_atom_set, rho_atom_set
141 :
142 1794 : CALL timeset(routineN, handle)
143 :
144 1794 : CPASSERT(ASSOCIATED(p_env%kpp1))
145 1794 : CPASSERT(ASSOCIATED(p_env%kpp1_env))
146 1794 : CPASSERT(ASSOCIATED(rho1))
147 :
148 1794 : my_calc_forces = .FALSE.
149 1794 : IF (PRESENT(calc_forces)) my_calc_forces = calc_forces
150 :
151 : CALL get_qs_env(qs_env, &
152 : pw_env=pw_env, &
153 : dft_control=dft_control, &
154 : para_env=para_env, &
155 1794 : rho=rho)
156 :
157 1794 : IF (lrigpw) THEN
158 : CALL get_qs_env(qs_env, &
159 : lri_env=lri_env, &
160 : lri_density=lri_density, &
161 0 : atomic_kind_set=atomic_kind_set)
162 : END IF
163 :
164 1794 : gapw = dft_control%qs_control%gapw
165 1794 : gapw_xc = dft_control%qs_control%gapw_xc
166 1794 : IF (gapw_xc) THEN
167 0 : CPASSERT(ASSOCIATED(rho1_xc))
168 : END IF
169 1794 : do_onecenter = gapw .OR. gapw_xc
170 1794 : nspins = dft_control%nspins
171 :
172 1794 : kpp1_env => p_env%kpp1_env
173 1794 : CALL kpp1_check_i_alloc(kpp1_env, qs_env, xc_section)
174 :
175 1794 : CALL qs_rho_get(rho, rho_ao=rho_ao)
176 1794 : CALL qs_rho_get(rho1, rho_g=rho1_g)
177 :
178 : ! gets the tmp grids
179 1794 : CPASSERT(ASSOCIATED(pw_env))
180 : CALL pw_env_get(pw_env, auxbas_pw_pool=auxbas_pw_pool, &
181 1794 : poisson_env=poisson_env)
182 1794 : CALL auxbas_pw_pool%create_pw(v_hartree_rspace)
183 :
184 1794 : IF (gapw .OR. gapw_xc) THEN
185 0 : CALL prepare_gapw_den(qs_env, p_env%local_rho_set, do_rho0=(.NOT. gapw_xc))
186 : END IF
187 :
188 : ! *** calculate the hartree potential on the total density ***
189 1794 : CALL auxbas_pw_pool%create_pw(rho1_tot_gspace)
190 1794 : CALL pw_zero(rho1_tot_gspace)
191 4148 : DO ispin = 1, nspins
192 4148 : CALL pw_axpy(rho1_g(ispin), rho1_tot_gspace)
193 : END DO
194 1794 : IF (gapw) THEN
195 0 : CALL pw_axpy(p_env%local_rho_set%rho0_mpole%rho0_s_gs, rho1_tot_gspace)
196 0 : IF (ASSOCIATED(p_env%local_rho_set%rho0_mpole%rhoz_cneo_s_gs)) THEN
197 0 : CALL pw_axpy(p_env%local_rho_set%rho0_mpole%rhoz_cneo_s_gs, rho1_tot_gspace)
198 : END IF
199 : END IF
200 : !
201 : BLOCK
202 : TYPE(pw_c1d_gs_type) :: v_hartree_gspace
203 1794 : CALL auxbas_pw_pool%create_pw(v_hartree_gspace)
204 : CALL pw_poisson_solve(poisson_env, rho1_tot_gspace, &
205 : energy_hartree, &
206 1794 : v_hartree_gspace)
207 1794 : CALL pw_transfer(v_hartree_gspace, v_hartree_rspace)
208 1794 : CALL auxbas_pw_pool%give_back_pw(v_hartree_gspace)
209 : END BLOCK
210 : !
211 1794 : CALL pw_scale(v_hartree_rspace, v_hartree_rspace%pw_grid%dvol)
212 :
213 1794 : CALL auxbas_pw_pool%give_back_pw(rho1_tot_gspace)
214 :
215 : ! *** calculate the xc potential ***
216 1794 : NULLIFY (v_xc, v_xc_tau)
217 1794 : IF (gapw_xc) THEN
218 0 : CALL get_qs_env(qs_env, rho_xc=rho0_fxc)
219 0 : rho1_fxc => rho1_xc
220 : ELSE
221 1794 : CALL get_qs_env(qs_env, rho=rho0_fxc)
222 1794 : rho1_fxc => rho1
223 : END IF
224 :
225 1794 : NULLIFY (rho_atom_set, rho1_atom_set)
226 1794 : IF (do_onecenter) THEN
227 0 : CALL get_qs_env(qs_env, rho_atom_set=rho_atom_set)
228 0 : rho1_atom_set => p_env%local_rho_set%rho_atom_set
229 : END IF
230 : CALL qs_fxc_apply(qs_env, kpp1_env%deriv_set, kpp1_env%rho_set, &
231 : rho1_fxc, rho_atom_set, xc_section, &
232 : do_onecenter, v_xc, v_xc_tau, rho1_atom_set, &
233 1794 : compute_virial=calc_virial, virial_xc=virial)
234 : ! nl-vdW
235 1794 : CALL qs_rho_get(rho0_fxc, rho_r=rho0_r)
236 1794 : CALL qs_rho_get(rho1_fxc, rho_r=rho1_r)
237 : CALL qs_fxc_nlvdw_apply(qs_env, xc_section, qs_env%dispersion_env, &
238 1794 : rho0_r, rho1_r, v_xc)
239 :
240 4148 : DO ispin = 1, nspins
241 4148 : CALL pw_scale(v_xc(ispin), v_xc(ispin)%pw_grid%dvol)
242 : END DO
243 1794 : v_rspace_new => v_xc
244 1794 : IF (SIZE(v_xc) /= nspins) THEN
245 0 : CALL auxbas_pw_pool%give_back_pw(v_xc(2))
246 : END IF
247 1794 : NULLIFY (v_xc)
248 1794 : IF (ASSOCIATED(v_xc_tau)) THEN
249 640 : DO ispin = 1, nspins
250 640 : CALL pw_scale(v_xc_tau(ispin), v_xc_tau(ispin)%pw_grid%dvol)
251 : END DO
252 252 : IF (SIZE(v_xc_tau) /= nspins) THEN
253 0 : CALL auxbas_pw_pool%give_back_pw(v_xc_tau(2))
254 : END IF
255 : END IF
256 :
257 1794 : alpha = 1.0_dp
258 1794 : IF (nspins == 1) alpha = 2.0_dp
259 :
260 : !-------------------------------!
261 : ! Add both hartree and xc terms !
262 : !-------------------------------!
263 4148 : DO ispin = 1, nspins
264 2354 : CALL dbcsr_set(kpp1_env%v_ao(ispin)%matrix, 0.0_dp)
265 :
266 2354 : IF (gapw_xc) THEN
267 : ! XC and Hartree are integrated separatedly
268 : ! XC uses the soft basis set only
269 : CALL integrate_v_rspace(v_rspace=v_rspace_new(ispin), &
270 : pmat=rho_ao(ispin), &
271 : hmat=kpp1_env%v_ao(ispin), &
272 : qs_env=qs_env, &
273 0 : calculate_forces=my_calc_forces, gapw=gapw_xc)
274 0 : IF (ASSOCIATED(v_xc_tau)) THEN
275 : CALL integrate_v_rspace(v_rspace=v_xc_tau(ispin), &
276 : pmat=rho_ao(ispin), &
277 : hmat=kpp1_env%v_ao(ispin), &
278 : qs_env=qs_env, &
279 : compute_tau=.TRUE., &
280 0 : calculate_forces=my_calc_forces, gapw=gapw_xc)
281 : END IF
282 : ! add Hartree for SINGLETS
283 0 : CALL pw_copy(v_hartree_rspace, v_rspace_new(1))
284 : CALL integrate_v_rspace(v_rspace=v_rspace_new(ispin), &
285 : pmat=rho_ao(ispin), &
286 : hmat=kpp1_env%v_ao(ispin), &
287 : qs_env=qs_env, &
288 0 : calculate_forces=my_calc_forces, gapw=gapw)
289 : ELSE
290 2354 : CALL pw_axpy(v_hartree_rspace, v_rspace_new(ispin))
291 2354 : IF (lrigpw) THEN
292 0 : IF (ASSOCIATED(v_xc_tau)) CPABORT("Meta-GGA functionals not supported with LRI!")
293 :
294 0 : lri_v_int => lri_density%lri_coefs(ispin)%lri_kinds
295 0 : CALL get_qs_env(qs_env, nkind=nkind)
296 0 : DO ikind = 1, nkind
297 0 : lri_v_int(ikind)%v_int = 0.0_dp
298 : END DO
299 : CALL integrate_v_rspace_one_center(v_rspace_new(ispin), qs_env, &
300 0 : lri_v_int, .FALSE., "LRI_AUX")
301 0 : DO ikind = 1, nkind
302 0 : CALL para_env%sum(lri_v_int(ikind)%v_int)
303 : END DO
304 0 : ALLOCATE (k1mat(1))
305 0 : k1mat(1)%matrix => kpp1_env%v_ao(ispin)%matrix
306 0 : IF (lri_env%exact_1c_terms) THEN
307 : CALL integrate_v_rspace_diagonal(v_rspace_new(ispin), k1mat(1)%matrix, &
308 0 : rho_ao(ispin)%matrix, qs_env, my_calc_forces, "ORB")
309 : END IF
310 0 : CALL calculate_lri_ks_matrix(lri_env, lri_v_int, k1mat, atomic_kind_set)
311 0 : DEALLOCATE (k1mat)
312 : ELSE
313 : CALL integrate_v_rspace(v_rspace=v_rspace_new(ispin), &
314 : pmat=rho_ao(ispin), &
315 : hmat=kpp1_env%v_ao(ispin), &
316 : qs_env=qs_env, &
317 2354 : calculate_forces=my_calc_forces, gapw=gapw)
318 2354 : IF (ASSOCIATED(v_xc_tau)) THEN
319 : CALL integrate_v_rspace(v_rspace=v_xc_tau(ispin), &
320 : pmat=rho_ao(ispin), &
321 : hmat=kpp1_env%v_ao(ispin), &
322 : qs_env=qs_env, &
323 : compute_tau=.TRUE., &
324 388 : calculate_forces=my_calc_forces, gapw=gapw)
325 : END IF
326 : END IF
327 : END IF
328 :
329 4148 : CALL dbcsr_add(p_env%kpp1(ispin)%matrix, kpp1_env%v_ao(ispin)%matrix, 1.0_dp, alpha)
330 : END DO
331 :
332 1794 : IF (gapw) THEN
333 : CALL Vh_1c_gg_integrals(qs_env, energy_hartree_1c, &
334 : p_env%hartree_local%ecoul_1c, &
335 : p_env%local_rho_set, &
336 0 : para_env, tddft=.TRUE., core_2nd=.TRUE.)
337 : CALL integrate_vhg0_rspace(qs_env, v_hartree_rspace, para_env, &
338 : calculate_forces=my_calc_forces, &
339 0 : local_rho_set=p_env%local_rho_set)
340 : ! *** Add single atom contributions to the KS matrix ***
341 : ! remap pointer
342 0 : ns = SIZE(p_env%kpp1)
343 0 : ksmat(1:ns, 1:1) => p_env%kpp1(1:ns)
344 0 : ns = SIZE(rho_ao)
345 0 : psmat(1:ns, 1:1) => rho_ao(1:ns)
346 : CALL update_ks_atom(qs_env, ksmat, psmat, forces=my_calc_forces, tddft=.TRUE., &
347 0 : rho_atom_external=p_env%local_rho_set%rho_atom_set)
348 1794 : ELSE IF (gapw_xc) THEN
349 0 : ns = SIZE(p_env%kpp1)
350 0 : ksmat(1:ns, 1:1) => p_env%kpp1(1:ns)
351 0 : ns = SIZE(rho_ao)
352 0 : psmat(1:ns, 1:1) => rho_ao(1:ns)
353 : CALL update_ks_atom(qs_env, ksmat, psmat, forces=my_calc_forces, tddft=.TRUE., &
354 0 : rho_atom_external=p_env%local_rho_set%rho_atom_set)
355 : END IF
356 :
357 1794 : CALL auxbas_pw_pool%give_back_pw(v_hartree_rspace)
358 4148 : DO ispin = 1, SIZE(v_rspace_new)
359 4148 : CALL auxbas_pw_pool%give_back_pw(v_rspace_new(ispin))
360 : END DO
361 1794 : DEALLOCATE (v_rspace_new)
362 1794 : IF (ASSOCIATED(v_xc_tau)) THEN
363 640 : DO ispin = 1, SIZE(v_xc_tau)
364 640 : CALL auxbas_pw_pool%give_back_pw(v_xc_tau(ispin))
365 : END DO
366 252 : DEALLOCATE (v_xc_tau)
367 : END IF
368 :
369 1794 : CALL timestop(handle)
370 :
371 3588 : END SUBROUTINE calc_kpp1
372 :
373 : ! **************************************************************************************************
374 : !> \brief checks that the intenal storage is allocated, and allocs it if needed
375 : !> \param kpp1_env the environment to check
376 : !> \param qs_env the qs environment this kpp1_env lives in
377 : !> \param xc_section ...
378 : !> \author Fawzi Mohamed
379 : !> \note
380 : !> private routine
381 : ! **************************************************************************************************
382 11562 : SUBROUTINE kpp1_check_i_alloc(kpp1_env, qs_env, xc_section)
383 :
384 : TYPE(qs_kpp1_env_type) :: kpp1_env
385 : TYPE(qs_environment_type), INTENT(IN), POINTER :: qs_env
386 : TYPE(section_vals_type), POINTER :: xc_section
387 :
388 : INTEGER :: ispin, nspins
389 : TYPE(admm_type), POINTER :: admm_env
390 11562 : TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_s
391 : TYPE(dft_control_type), POINTER :: dft_control
392 : TYPE(pw_env_type), POINTER :: pw_env
393 : TYPE(qs_rho_type), POINTER :: rho
394 : TYPE(section_vals_type), POINTER :: admm_xc_section
395 :
396 : ! ------------------------------------------------------------------
397 :
398 : CALL get_qs_env(qs_env, pw_env=pw_env, matrix_s=matrix_s, &
399 11562 : admm_env=admm_env, dft_control=dft_control)
400 :
401 11562 : nspins = dft_control%nspins
402 :
403 11562 : IF (.NOT. ASSOCIATED(kpp1_env%v_ao)) THEN
404 1504 : CALL dbcsr_allocate_matrix_set(kpp1_env%v_ao, nspins)
405 3242 : DO ispin = 1, nspins
406 1738 : ALLOCATE (kpp1_env%v_ao(ispin)%matrix)
407 : CALL dbcsr_copy(kpp1_env%v_ao(ispin)%matrix, matrix_s(1)%matrix, &
408 3242 : name="kpp1%v_ao-"//ADJUSTL(cp_to_string(ispin)))
409 : END DO
410 : END IF
411 :
412 11562 : IF (.NOT. ASSOCIATED(kpp1_env%deriv_set)) THEN
413 1504 : IF (dft_control%qs_control%gapw_xc) THEN
414 52 : CALL get_qs_env(qs_env, rho_xc=rho)
415 : ELSE
416 1452 : CALL get_qs_env(qs_env, rho=rho)
417 : END IF
418 34592 : ALLOCATE (kpp1_env%deriv_set, kpp1_env%rho_set)
419 : CALL qs_fxc_prep(qs_env, rho, &
420 : kpp1_env%rho_set, kpp1_env%deriv_set, &
421 1504 : xc_section, pw_env, is_triplet=.FALSE.)
422 : END IF
423 :
424 : ! ADMM Correction
425 11562 : IF (dft_control%do_admm) THEN
426 2370 : IF (admm_env%aux_exch_func /= do_admm_aux_exch_func_none) THEN
427 1404 : IF (.NOT. ASSOCIATED(kpp1_env%deriv_set_admm)) THEN
428 4462 : ALLOCATE (kpp1_env%deriv_set_admm, kpp1_env%rho_set_admm)
429 194 : admm_xc_section => admm_env%xc_section_aux
430 194 : CALL get_admm_env(admm_env, rho_aux_fit=rho)
431 : CALL qs_fxc_prep(qs_env, rho, kpp1_env%rho_set_admm, kpp1_env%deriv_set_admm, &
432 194 : admm_xc_section, pw_env, is_triplet=.FALSE.)
433 : END IF
434 : END IF
435 : END IF
436 :
437 11562 : END SUBROUTINE kpp1_check_i_alloc
438 : ! **************************************************************************************************
439 : !> \brief ...
440 : !> \param rho1 ...
441 : !> \param rho1_tot_gspace ...
442 : !> \param out_unit ...
443 : ! **************************************************************************************************
444 0 : SUBROUTINE print_densities(rho1, rho1_tot_gspace, out_unit)
445 :
446 : TYPE(qs_rho_type), POINTER :: rho1
447 : TYPE(pw_c1d_gs_type), INTENT(IN) :: rho1_tot_gspace
448 : INTEGER :: out_unit
449 :
450 : REAL(KIND=dp) :: total_rho_gspace
451 0 : REAL(KIND=dp), DIMENSION(:), POINTER :: tot_rho1_r
452 :
453 0 : NULLIFY (tot_rho1_r)
454 :
455 0 : total_rho_gspace = pw_integrate_function(rho1_tot_gspace, isign=-1)
456 0 : IF (out_unit > 0) THEN
457 0 : CALL qs_rho_get(rho1, tot_rho_r=tot_rho1_r)
458 : WRITE (UNIT=out_unit, FMT="(T3,A,T60,F20.10)") &
459 0 : "KPP1 total charge density (r-space):", &
460 0 : accurate_sum(tot_rho1_r), &
461 0 : "KPP1 total charge density (g-space):", &
462 0 : total_rho_gspace
463 : END IF
464 :
465 0 : END SUBROUTINE print_densities
466 :
467 : END MODULE qs_kpp1_env_methods
|
