(git:15c1bfc)
Loading...
Searching...
No Matches
qs_kpp1_env_methods.F
Go to the documentation of this file.
1!--------------------------------------------------------------------------------------------------!
2! CP2K: A general program to perform molecular dynamics simulations !
3! Copyright 2000-2025 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! **************************************************************************************************
16 USE admm_types, ONLY: admm_type,&
20 USE cp_dbcsr_api, ONLY: dbcsr_add,&
38 USE kahan_sum, ONLY: accurate_sum
39 USE kinds, ONLY: dp
45 USE pw_env_types, ONLY: pw_env_get,&
47 USE pw_methods, ONLY: pw_axpy,&
48 pw_copy,&
50 pw_scale,&
55 USE pw_types, ONLY: pw_c1d_gs_type,&
60 USE qs_integrate_potential, ONLY: integrate_v_rspace,&
61 integrate_v_rspace_diagonal,&
62 integrate_v_rspace_one_center
64 USE qs_ks_atom, ONLY: update_ks_atom
68 USE qs_rho_types, ONLY: qs_rho_get,&
71 USE xc, ONLY: xc_calc_2nd_deriv,&
75#include "./base/base_uses.f90"
76
77 IMPLICIT NONE
78
79 PRIVATE
80
81 LOGICAL, PRIVATE, PARAMETER :: debug_this_module = .true.
82 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'qs_kpp1_env_methods'
83
84 PUBLIC :: kpp1_create, &
87
88CONTAINS
89
90! **************************************************************************************************
91!> \brief allocates and initializes a kpp1_env
92!> \param kpp1_env the environment to initialize
93!> \par History
94!> 07.2002 created [fawzi]
95!> \author Fawzi Mohamed
96! **************************************************************************************************
97 SUBROUTINE kpp1_create(kpp1_env)
98 TYPE(qs_kpp1_env_type) :: kpp1_env
99
100 NULLIFY (kpp1_env%v_ao, kpp1_env%rho_set, kpp1_env%deriv_set, &
101 kpp1_env%rho_set_admm, kpp1_env%deriv_set_admm)
102 END SUBROUTINE kpp1_create
103
104! **************************************************************************************************
105!> \brief ...
106!> \param rho1_xc ...
107!> \param rho1 ...
108!> \param xc_section ...
109!> \param lrigpw ...
110!> \param do_triplet ...
111!> \param qs_env ...
112!> \param p_env ...
113!> \param calc_forces ...
114!> \param calc_virial ...
115!> \param virial ...
116! **************************************************************************************************
117 SUBROUTINE calc_kpp1(rho1_xc, rho1, xc_section, lrigpw, do_triplet, qs_env, p_env, &
118 calc_forces, calc_virial, virial)
119
120 TYPE(qs_rho_type), POINTER :: rho1_xc, rho1
121 TYPE(section_vals_type), POINTER :: xc_section
122 LOGICAL, INTENT(IN) :: lrigpw, do_triplet
123 TYPE(qs_environment_type), POINTER :: qs_env
124 TYPE(qs_p_env_type) :: p_env
125 LOGICAL, INTENT(IN), OPTIONAL :: calc_forces, calc_virial
126 REAL(kind=dp), DIMENSION(3, 3), INTENT(INOUT), &
127 OPTIONAL :: virial
128
129 CHARACTER(len=*), PARAMETER :: routinen = 'calc_kpp1'
130
131 INTEGER :: handle, ikind, ispin, nkind, ns, nspins, &
132 output_unit
133 LOGICAL :: gapw, gapw_xc, lsd, my_calc_forces
134 REAL(kind=dp) :: alpha, energy_hartree, energy_hartree_1c
135 TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
136 TYPE(cp_logger_type), POINTER :: logger
137 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: k1mat, rho_ao
138 TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: ksmat, psmat
139 TYPE(lri_density_type), POINTER :: lri_density
140 TYPE(lri_environment_type), POINTER :: lri_env
141 TYPE(lri_kind_type), DIMENSION(:), POINTER :: lri_v_int
142 TYPE(mp_para_env_type), POINTER :: para_env
143 TYPE(pw_c1d_gs_type) :: rho1_tot_gspace
144 TYPE(pw_c1d_gs_type), DIMENSION(:), POINTER :: rho1_g, rho1_g_pw
145 TYPE(pw_env_type), POINTER :: pw_env
146 TYPE(pw_poisson_type), POINTER :: poisson_env
147 TYPE(pw_pool_type), POINTER :: auxbas_pw_pool
148 TYPE(pw_r3d_rs_type) :: v_hartree_rspace
149 TYPE(pw_r3d_rs_type), DIMENSION(:), POINTER :: rho1_r, rho1_r_pw, tau1_r, tau1_r_pw, &
150 v_rspace_new, v_xc, v_xc_tau
151 TYPE(qs_rho_type), POINTER :: rho
152 TYPE(rho_atom_type), DIMENSION(:), POINTER :: rho1_atom_set, rho_atom_set
153 TYPE(section_vals_type), POINTER :: input, scf_section
154
155 CALL timeset(routinen, handle)
156
157 NULLIFY (v_xc, rho1_g, pw_env, rho1_g_pw, tau1_r_pw)
158 logger => cp_get_default_logger()
159
160 cpassert(ASSOCIATED(p_env%kpp1))
161 cpassert(ASSOCIATED(p_env%kpp1_env))
162 cpassert(ASSOCIATED(rho1))
163
164 nspins = SIZE(p_env%kpp1)
165 lsd = (nspins == 2)
166
167 my_calc_forces = .false.
168 IF (PRESENT(calc_forces)) my_calc_forces = calc_forces
169
170 CALL get_qs_env(qs_env, &
171 pw_env=pw_env, &
172 input=input, &
173 para_env=para_env, &
174 rho=rho)
175
176 cpassert(ASSOCIATED(rho1))
177
178 IF (lrigpw) THEN
179 CALL get_qs_env(qs_env, &
180 lri_env=lri_env, &
181 lri_density=lri_density, &
182 atomic_kind_set=atomic_kind_set)
183 END IF
184
185 gapw = (section_get_ival(input, "DFT%QS%METHOD") == do_method_gapw)
186 gapw_xc = (section_get_ival(input, "DFT%QS%METHOD") == do_method_gapw_xc)
187 IF (gapw_xc) THEN
188 cpassert(ASSOCIATED(rho1_xc))
189 END IF
190
191 CALL kpp1_check_i_alloc(p_env%kpp1_env, qs_env, do_triplet)
192
193 CALL qs_rho_get(rho, rho_ao=rho_ao)
194 CALL qs_rho_get(rho1, rho_g=rho1_g)
195
196 ! gets the tmp grids
197 cpassert(ASSOCIATED(pw_env))
198 CALL pw_env_get(pw_env, auxbas_pw_pool=auxbas_pw_pool, &
199 poisson_env=poisson_env)
200 CALL auxbas_pw_pool%create_pw(v_hartree_rspace)
201
202 IF (gapw .OR. gapw_xc) &
203 CALL prepare_gapw_den(qs_env, p_env%local_rho_set, do_rho0=(.NOT. gapw_xc))
204
205 ! *** calculate the hartree potential on the total density ***
206 CALL auxbas_pw_pool%create_pw(rho1_tot_gspace)
207
208 CALL pw_copy(rho1_g(1), rho1_tot_gspace)
209 DO ispin = 2, nspins
210 CALL pw_axpy(rho1_g(ispin), rho1_tot_gspace)
211 END DO
212 IF (gapw) THEN
213 CALL pw_axpy(p_env%local_rho_set%rho0_mpole%rho0_s_gs, rho1_tot_gspace)
214 IF (ASSOCIATED(p_env%local_rho_set%rho0_mpole%rhoz_cneo_s_gs)) THEN
215 CALL pw_axpy(p_env%local_rho_set%rho0_mpole%rhoz_cneo_s_gs, rho1_tot_gspace)
216 END IF
217 END IF
218
219 scf_section => section_vals_get_subs_vals(input, "DFT%SCF")
220 IF (cp_print_key_should_output(logger%iter_info, scf_section, "PRINT%TOTAL_DENSITIES") &
221 /= 0) THEN
222 output_unit = cp_print_key_unit_nr(logger, scf_section, "PRINT%TOTAL_DENSITIES", &
223 extension=".scfLog")
224 CALL print_densities(rho1, rho1_tot_gspace, output_unit)
225 CALL cp_print_key_finished_output(output_unit, logger, scf_section, &
226 "PRINT%TOTAL_DENSITIES")
227 END IF
228
229 IF (.NOT. (nspins == 1 .AND. do_triplet)) THEN
230 block
231 TYPE(pw_c1d_gs_type) :: v_hartree_gspace
232 CALL auxbas_pw_pool%create_pw(v_hartree_gspace)
233 CALL pw_poisson_solve(poisson_env, rho1_tot_gspace, &
234 energy_hartree, &
235 v_hartree_gspace)
236 CALL pw_transfer(v_hartree_gspace, v_hartree_rspace)
237 CALL auxbas_pw_pool%give_back_pw(v_hartree_gspace)
238 END block
239 CALL pw_scale(v_hartree_rspace, v_hartree_rspace%pw_grid%dvol)
240 END IF
241
242 CALL auxbas_pw_pool%give_back_pw(rho1_tot_gspace)
243
244 ! *** calculate the xc potential ***
245 IF (gapw_xc) THEN
246 CALL qs_rho_get(rho1_xc, rho_r=rho1_r, tau_r=tau1_r)
247 ELSE
248 CALL qs_rho_get(rho1, rho_r=rho1_r, tau_r=tau1_r)
249 END IF
250
251 IF (nspins == 1 .AND. do_triplet) THEN
252
253 lsd = .true.
254 ALLOCATE (rho1_r_pw(2))
255 DO ispin = 1, 2
256 CALL rho1_r_pw(ispin)%create(rho1_r(1)%pw_grid)
257 CALL pw_transfer(rho1_r(1), rho1_r_pw(ispin))
258 END DO
259
260 IF (ASSOCIATED(tau1_r)) THEN
261 ALLOCATE (tau1_r_pw(2))
262 DO ispin = 1, 2
263 CALL tau1_r_pw(ispin)%create(tau1_r(1)%pw_grid)
264 CALL pw_transfer(tau1_r(1), tau1_r_pw(ispin))
265 END DO
266 END IF
267
268 ELSE
269
270 rho1_r_pw => rho1_r
271
272 tau1_r_pw => tau1_r
273
274 END IF
275
276 CALL xc_calc_2nd_deriv(v_xc, v_xc_tau, p_env%kpp1_env%deriv_set, p_env%kpp1_env%rho_set, &
277 rho1_r_pw, rho1_g_pw, tau1_r_pw, auxbas_pw_pool, xc_section, .false., &
278 do_excitations=.true., do_triplet=do_triplet, &
279 compute_virial=calc_virial, virial_xc=virial)
280
281 DO ispin = 1, nspins
282 CALL pw_scale(v_xc(ispin), v_xc(ispin)%pw_grid%dvol)
283 END DO
284 v_rspace_new => v_xc
285 IF (SIZE(v_xc) /= nspins) THEN
286 CALL auxbas_pw_pool%give_back_pw(v_xc(2))
287 END IF
288 NULLIFY (v_xc)
289 IF (ASSOCIATED(v_xc_tau)) THEN
290 DO ispin = 1, nspins
291 CALL pw_scale(v_xc_tau(ispin), v_xc_tau(ispin)%pw_grid%dvol)
292 END DO
293 IF (SIZE(v_xc_tau) /= nspins) THEN
294 CALL auxbas_pw_pool%give_back_pw(v_xc_tau(2))
295 END IF
296 END IF
297
298 IF (gapw .OR. gapw_xc) THEN
299 CALL get_qs_env(qs_env, rho_atom_set=rho_atom_set)
300 rho1_atom_set => p_env%local_rho_set%rho_atom_set
301 CALL calculate_xc_2nd_deriv_atom(rho_atom_set, rho1_atom_set, qs_env, xc_section, para_env, &
302 do_triplet=do_triplet)
303 END IF
304
305 IF (nspins == 1 .AND. do_triplet) THEN
306 DO ispin = 1, SIZE(rho1_r_pw)
307 CALL rho1_r_pw(ispin)%release()
308 END DO
309 DEALLOCATE (rho1_r_pw)
310 IF (ASSOCIATED(tau1_r_pw)) THEN
311 DO ispin = 1, SIZE(tau1_r_pw)
312 CALL tau1_r_pw(ispin)%release()
313 END DO
314 DEALLOCATE (tau1_r_pw)
315 END IF
316 END IF
317
318 alpha = 1.0_dp
319 IF (nspins == 1) alpha = 2.0_dp
320
321 !-------------------------------!
322 ! Add both hartree and xc terms !
323 !-------------------------------!
324 DO ispin = 1, nspins
325 CALL dbcsr_set(p_env%kpp1_env%v_ao(ispin)%matrix, 0.0_dp)
326
327 ! XC and Hartree are integrated separatedly
328 ! XC uses the soft basis set only
329 IF (gapw_xc) THEN
330
331 IF (nspins == 1) THEN
332 CALL integrate_v_rspace(v_rspace=v_rspace_new(ispin), &
333 pmat=rho_ao(ispin), &
334 hmat=p_env%kpp1_env%v_ao(ispin), &
335 qs_env=qs_env, &
336 calculate_forces=my_calc_forces, gapw=gapw_xc)
337
338 IF (ASSOCIATED(v_xc_tau)) THEN
339 CALL integrate_v_rspace(v_rspace=v_xc_tau(ispin), &
340 pmat=rho_ao(ispin), &
341 hmat=p_env%kpp1_env%v_ao(ispin), &
342 qs_env=qs_env, &
343 compute_tau=.true., &
344 calculate_forces=my_calc_forces, gapw=gapw_xc)
345 END IF
346
347 ! add hartree only for SINGLETS
348 IF (.NOT. do_triplet) THEN
349 CALL pw_copy(v_hartree_rspace, v_rspace_new(1))
350
351 CALL integrate_v_rspace(v_rspace=v_rspace_new(ispin), &
352 pmat=rho_ao(ispin), &
353 hmat=p_env%kpp1_env%v_ao(ispin), &
354 qs_env=qs_env, &
355 calculate_forces=my_calc_forces, gapw=gapw)
356 END IF
357 ELSE
358 CALL integrate_v_rspace(v_rspace=v_rspace_new(ispin), &
359 pmat=rho_ao(ispin), &
360 hmat=p_env%kpp1_env%v_ao(ispin), &
361 qs_env=qs_env, &
362 calculate_forces=my_calc_forces, gapw=gapw_xc)
363
364 IF (ASSOCIATED(v_xc_tau)) THEN
365 CALL integrate_v_rspace(v_rspace=v_xc_tau(ispin), &
366 pmat=rho_ao(ispin), &
367 hmat=p_env%kpp1_env%v_ao(ispin), &
368 qs_env=qs_env, &
369 compute_tau=.true., &
370 calculate_forces=my_calc_forces, gapw=gapw_xc)
371 END IF
372
373 CALL pw_copy(v_hartree_rspace, v_rspace_new(ispin))
374 CALL integrate_v_rspace(v_rspace=v_rspace_new(ispin), &
375 pmat=rho_ao(ispin), &
376 hmat=p_env%kpp1_env%v_ao(ispin), &
377 qs_env=qs_env, &
378 calculate_forces=my_calc_forces, gapw=gapw)
379 END IF
380
381 ELSE
382
383 IF (nspins == 1) THEN
384
385 ! add hartree only for SINGLETS
386 IF (.NOT. do_triplet) THEN
387 CALL pw_axpy(v_hartree_rspace, v_rspace_new(1))
388 END IF
389 ELSE
390 CALL pw_axpy(v_hartree_rspace, v_rspace_new(ispin))
391 END IF
392
393 IF (lrigpw) THEN
394 IF (ASSOCIATED(v_xc_tau)) cpabort("Meta-GGA functionals not supported with LRI!")
395
396 lri_v_int => lri_density%lri_coefs(ispin)%lri_kinds
397 CALL get_qs_env(qs_env, nkind=nkind)
398 DO ikind = 1, nkind
399 lri_v_int(ikind)%v_int = 0.0_dp
400 END DO
401 CALL integrate_v_rspace_one_center(v_rspace_new(ispin), qs_env, &
402 lri_v_int, .false., "LRI_AUX")
403 DO ikind = 1, nkind
404 CALL para_env%sum(lri_v_int(ikind)%v_int)
405 END DO
406 ALLOCATE (k1mat(1))
407 k1mat(1)%matrix => p_env%kpp1_env%v_ao(ispin)%matrix
408 IF (lri_env%exact_1c_terms) THEN
409 CALL integrate_v_rspace_diagonal(v_rspace_new(ispin), k1mat(1)%matrix, &
410 rho_ao(ispin)%matrix, qs_env, my_calc_forces, "ORB")
411 END IF
412 CALL calculate_lri_ks_matrix(lri_env, lri_v_int, k1mat, atomic_kind_set)
413 DEALLOCATE (k1mat)
414 ELSE
415 CALL integrate_v_rspace(v_rspace=v_rspace_new(ispin), &
416 pmat=rho_ao(ispin), &
417 hmat=p_env%kpp1_env%v_ao(ispin), &
418 qs_env=qs_env, &
419 calculate_forces=my_calc_forces, gapw=gapw)
420
421 IF (ASSOCIATED(v_xc_tau)) THEN
422 CALL integrate_v_rspace(v_rspace=v_xc_tau(ispin), &
423 pmat=rho_ao(ispin), &
424 hmat=p_env%kpp1_env%v_ao(ispin), &
425 qs_env=qs_env, &
426 compute_tau=.true., &
427 calculate_forces=my_calc_forces, gapw=gapw)
428 END IF
429 END IF
430 END IF
431
432 CALL dbcsr_add(p_env%kpp1(ispin)%matrix, p_env%kpp1_env%v_ao(ispin)%matrix, 1.0_dp, alpha)
433 END DO
434
435 IF (gapw) THEN
436 IF (.NOT. (nspins == 1 .AND. do_triplet)) THEN
437 CALL vh_1c_gg_integrals(qs_env, energy_hartree_1c, &
438 p_env%hartree_local%ecoul_1c, &
439 p_env%local_rho_set, &
440 para_env, tddft=.true., core_2nd=.true.)
441 CALL integrate_vhg0_rspace(qs_env, v_hartree_rspace, para_env, &
442 calculate_forces=my_calc_forces, &
443 local_rho_set=p_env%local_rho_set)
444 END IF
445 ! *** Add single atom contributions to the KS matrix ***
446 ! remap pointer
447 ns = SIZE(p_env%kpp1)
448 ksmat(1:ns, 1:1) => p_env%kpp1(1:ns)
449 ns = SIZE(rho_ao)
450 psmat(1:ns, 1:1) => rho_ao(1:ns)
451 CALL update_ks_atom(qs_env, ksmat, psmat, forces=my_calc_forces, tddft=.true., &
452 rho_atom_external=p_env%local_rho_set%rho_atom_set)
453 ELSEIF (gapw_xc) THEN
454 ns = SIZE(p_env%kpp1)
455 ksmat(1:ns, 1:1) => p_env%kpp1(1:ns)
456 ns = SIZE(rho_ao)
457 psmat(1:ns, 1:1) => rho_ao(1:ns)
458 CALL update_ks_atom(qs_env, ksmat, psmat, forces=my_calc_forces, tddft=.true., &
459 rho_atom_external=p_env%local_rho_set%rho_atom_set)
460 END IF
461
462 CALL auxbas_pw_pool%give_back_pw(v_hartree_rspace)
463 DO ispin = 1, SIZE(v_rspace_new)
464 CALL auxbas_pw_pool%give_back_pw(v_rspace_new(ispin))
465 END DO
466 DEALLOCATE (v_rspace_new)
467 IF (ASSOCIATED(v_xc_tau)) THEN
468 DO ispin = 1, SIZE(v_xc_tau)
469 CALL auxbas_pw_pool%give_back_pw(v_xc_tau(ispin))
470 END DO
471 DEALLOCATE (v_xc_tau)
472 END IF
473
474 CALL timestop(handle)
475 END SUBROUTINE calc_kpp1
476
477! **************************************************************************************************
478!> \brief checks that the intenal storage is allocated, and allocs it if needed
479!> \param kpp1_env the environment to check
480!> \param qs_env the qs environment this kpp1_env lives in
481!> \param do_triplet ...
482!> \author Fawzi Mohamed
483!> \note
484!> private routine
485! **************************************************************************************************
486 SUBROUTINE kpp1_check_i_alloc(kpp1_env, qs_env, do_triplet)
487
488 TYPE(qs_kpp1_env_type) :: kpp1_env
489 TYPE(qs_environment_type), INTENT(IN), POINTER :: qs_env
490 LOGICAL, INTENT(IN) :: do_triplet
491
492 INTEGER :: ispin, nspins
493 TYPE(admm_type), POINTER :: admm_env
494 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_s
495 TYPE(dft_control_type), POINTER :: dft_control
496 TYPE(pw_env_type), POINTER :: pw_env
497 TYPE(pw_pool_type), POINTER :: auxbas_pw_pool
498 TYPE(pw_r3d_rs_type), DIMENSION(:), POINTER :: my_rho_r, my_tau_r, rho_r, tau_r
499 TYPE(qs_rho_type), POINTER :: rho
500 TYPE(section_vals_type), POINTER :: admm_xc_section, input, xc_section
501
502! ------------------------------------------------------------------
503
504 NULLIFY (pw_env, auxbas_pw_pool, matrix_s, rho, rho_r, admm_env, dft_control, my_rho_r, my_tau_r)
505
506 CALL get_qs_env(qs_env, pw_env=pw_env, &
507 matrix_s=matrix_s, rho=rho, input=input, &
508 admm_env=admm_env, dft_control=dft_control)
509
510 CALL qs_rho_get(rho, rho_r=rho_r, tau_r=tau_r)
511 nspins = SIZE(rho_r)
512
513 CALL pw_env_get(pw_env, auxbas_pw_pool=auxbas_pw_pool)
514
515 IF (.NOT. ASSOCIATED(kpp1_env%v_ao)) THEN
516 CALL dbcsr_allocate_matrix_set(kpp1_env%v_ao, nspins)
517 DO ispin = 1, nspins
518 ALLOCATE (kpp1_env%v_ao(ispin)%matrix)
519 CALL dbcsr_copy(kpp1_env%v_ao(ispin)%matrix, matrix_s(1)%matrix, &
520 name="kpp1%v_ao-"//adjustl(cp_to_string(ispin)))
521 END DO
522 END IF
523
524 IF (.NOT. ASSOCIATED(kpp1_env%deriv_set)) THEN
525
526 IF (nspins == 1 .AND. do_triplet) THEN
527 ALLOCATE (my_rho_r(2))
528 DO ispin = 1, 2
529 CALL auxbas_pw_pool%create_pw(my_rho_r(ispin))
530 CALL pw_axpy(rho_r(1), my_rho_r(ispin), 0.5_dp, 0.0_dp)
531 END DO
532 IF (dft_control%use_kinetic_energy_density) THEN
533 ALLOCATE (my_tau_r(2))
534 DO ispin = 1, 2
535 CALL auxbas_pw_pool%create_pw(my_tau_r(ispin))
536 CALL pw_axpy(tau_r(1), my_tau_r(ispin), 0.5_dp, 0.0_dp)
537 END DO
538 END IF
539 ELSE
540 my_rho_r => rho_r
541 IF (dft_control%use_kinetic_energy_density) THEN
542 my_tau_r => tau_r
543 END IF
544 END IF
545
546 IF (dft_control%do_admm) THEN
547 xc_section => admm_env%xc_section_primary
548 ELSE
549 xc_section => section_vals_get_subs_vals(input, "DFT%XC")
550 END IF
551
552 ALLOCATE (kpp1_env%deriv_set, kpp1_env%rho_set)
553 CALL xc_prep_2nd_deriv(kpp1_env%deriv_set, kpp1_env%rho_set, &
554 my_rho_r, auxbas_pw_pool, &
555 xc_section=xc_section, tau_r=my_tau_r)
556
557 IF (nspins == 1 .AND. do_triplet) THEN
558 DO ispin = 1, SIZE(my_rho_r)
559 CALL my_rho_r(ispin)%release()
560 END DO
561 DEALLOCATE (my_rho_r)
562 IF (ASSOCIATED(my_tau_r)) THEN
563 DO ispin = 1, SIZE(my_tau_r)
564 CALL my_tau_r(ispin)%release()
565 END DO
566 DEALLOCATE (my_tau_r)
567 END IF
568 END IF
569 END IF
570
571 ! ADMM Correction
572 IF (dft_control%do_admm) THEN
573 IF (admm_env%aux_exch_func /= do_admm_aux_exch_func_none) THEN
574 IF (.NOT. ASSOCIATED(kpp1_env%deriv_set_admm)) THEN
575 cpassert(.NOT. do_triplet)
576 admm_xc_section => admm_env%xc_section_aux
577 CALL get_admm_env(qs_env%admm_env, rho_aux_fit=rho)
578 CALL qs_rho_get(rho, rho_r=rho_r)
579 ALLOCATE (kpp1_env%deriv_set_admm, kpp1_env%rho_set_admm)
580 CALL xc_prep_2nd_deriv(kpp1_env%deriv_set_admm, kpp1_env%rho_set_admm, &
581 rho_r, auxbas_pw_pool, &
582 xc_section=admm_xc_section)
583 END IF
584 END IF
585 END IF
586
587 END SUBROUTINE kpp1_check_i_alloc
588
589! **************************************************************************************************
590!> \brief function to advise of changes either in the grids
591!> \param kpp1_env the kpp1_env
592!> \par History
593!> 11.2002 created [fawzi]
594!> \author Fawzi Mohamed
595! **************************************************************************************************
596 SUBROUTINE kpp1_did_change(kpp1_env)
597 TYPE(qs_kpp1_env_type) :: kpp1_env
598
599 IF (ASSOCIATED(kpp1_env%deriv_set)) THEN
600 CALL xc_dset_release(kpp1_env%deriv_set)
601 DEALLOCATE (kpp1_env%deriv_set)
602 NULLIFY (kpp1_env%deriv_set)
603 END IF
604 IF (ASSOCIATED(kpp1_env%rho_set)) THEN
605 CALL xc_rho_set_release(kpp1_env%rho_set)
606 DEALLOCATE (kpp1_env%rho_set)
607 END IF
608
609 END SUBROUTINE kpp1_did_change
610
611! **************************************************************************************************
612!> \brief ...
613!> \param rho1 ...
614!> \param rho1_tot_gspace ...
615!> \param out_unit ...
616! **************************************************************************************************
617 SUBROUTINE print_densities(rho1, rho1_tot_gspace, out_unit)
618
619 TYPE(qs_rho_type), POINTER :: rho1
620 TYPE(pw_c1d_gs_type), INTENT(IN) :: rho1_tot_gspace
621 INTEGER :: out_unit
622
623 REAL(kind=dp) :: total_rho_gspace
624 REAL(kind=dp), DIMENSION(:), POINTER :: tot_rho1_r
625
626 NULLIFY (tot_rho1_r)
627
628 total_rho_gspace = pw_integrate_function(rho1_tot_gspace, isign=-1)
629 IF (out_unit > 0) THEN
630 CALL qs_rho_get(rho1, tot_rho_r=tot_rho1_r)
631 WRITE (unit=out_unit, fmt="(T3,A,T60,F20.10)") &
632 "KPP1 total charge density (r-space):", &
633 accurate_sum(tot_rho1_r), &
634 "KPP1 total charge density (g-space):", &
635 total_rho_gspace
636 END IF
637
638 END SUBROUTINE print_densities
639
640END MODULE qs_kpp1_env_methods
Types and set/get functions for auxiliary density matrix methods.
Definition admm_types.F:15
subroutine, public get_admm_env(admm_env, mo_derivs_aux_fit, mos_aux_fit, sab_aux_fit, sab_aux_fit_asymm, sab_aux_fit_vs_orb, matrix_s_aux_fit, matrix_s_aux_fit_kp, matrix_s_aux_fit_vs_orb, matrix_s_aux_fit_vs_orb_kp, task_list_aux_fit, matrix_ks_aux_fit, matrix_ks_aux_fit_kp, matrix_ks_aux_fit_im, matrix_ks_aux_fit_dft, matrix_ks_aux_fit_hfx, matrix_ks_aux_fit_dft_kp, matrix_ks_aux_fit_hfx_kp, rho_aux_fit, rho_aux_fit_buffer, admm_dm)
Get routine for the ADMM env.
Definition admm_types.F:593
Define the atomic kind types and their sub types.
Defines control structures, which contain the parameters and the settings for the DFT-based calculati...
subroutine, public dbcsr_copy(matrix_b, matrix_a, name, keep_sparsity, keep_imaginary)
...
subroutine, public dbcsr_set(matrix, alpha)
...
subroutine, public dbcsr_add(matrix_a, matrix_b, alpha_scalar, beta_scalar)
...
DBCSR operations in CP2K.
various routines to log and control the output. The idea is that decisions about where to log should ...
type(cp_logger_type) function, pointer, public cp_get_default_logger()
returns the default logger
routines to handle the output, The idea is to remove the decision of wheter to output and what to out...
integer function, public cp_print_key_unit_nr(logger, basis_section, print_key_path, extension, middle_name, local, log_filename, ignore_should_output, file_form, file_position, file_action, file_status, do_backup, on_file, is_new_file, mpi_io, fout)
...
subroutine, public cp_print_key_finished_output(unit_nr, logger, basis_section, print_key_path, local, ignore_should_output, on_file, mpi_io)
should be called after you finish working with a unit obtained with cp_print_key_unit_nr,...
integer function, public cp_print_key_should_output(iteration_info, basis_section, print_key_path, used_print_key, first_time)
returns what should be done with the given property if btest(res,cp_p_store) then the property should...
subroutine, public vh_1c_gg_integrals(qs_env, energy_hartree_1c, ecoul_1c, local_rho_set, para_env, tddft, local_rho_set_2nd, core_2nd)
Calculates one center GAPW Hartree energies and matrix elements Hartree potentials are input Takes po...
collects all constants needed in input so that they can be used without circular dependencies
integer, parameter, public do_admm_aux_exch_func_none
integer, parameter, public do_method_gapw
integer, parameter, public do_method_gapw_xc
objects that represent the structure of input sections and the data contained in an input section
integer function, public section_get_ival(section_vals, keyword_name)
...
recursive type(section_vals_type) function, pointer, public section_vals_get_subs_vals(section_vals, subsection_name, i_rep_section, can_return_null)
returns the values of the requested subsection
sums arrays of real/complex numbers with much reduced round-off as compared to a naive implementation...
Definition kahan_sum.F:29
Defines the basic variable types.
Definition kinds.F:23
integer, parameter, public dp
Definition kinds.F:34
contains the types and subroutines for dealing with the lri_env lri : local resolution of the identit...
routines that build the Kohn-Sham matrix for the LRIGPW and xc parts
subroutine, public calculate_lri_ks_matrix(lri_env, lri_v_int, h_matrix, atomic_kind_set, cell_to_index)
update of LRIGPW KS matrix
Interface to the message passing library MPI.
container for various plainwaves related things
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
functions related to the poisson solver on regular grids
Manages a pool of grids (to be used for example as tmp objects), but can also be used to instantiate ...
subroutine, public get_qs_env(qs_env, atomic_kind_set, qs_kind_set, cell, super_cell, cell_ref, use_ref_cell, kpoints, dft_control, mos, sab_orb, sab_all, qmmm, qmmm_periodic, sac_ae, sac_ppl, sac_lri, sap_ppnl, sab_vdw, sab_scp, sap_oce, sab_lrc, sab_se, sab_xtbe, sab_tbe, sab_core, sab_xb, sab_xtb_pp, sab_xtb_nonbond, sab_almo, sab_kp, sab_kp_nosym, sab_cneo, particle_set, energy, force, matrix_h, matrix_h_im, matrix_ks, matrix_ks_im, matrix_vxc, run_rtp, rtp, matrix_h_kp, matrix_h_im_kp, matrix_ks_kp, matrix_ks_im_kp, matrix_vxc_kp, kinetic_kp, matrix_s_kp, matrix_w_kp, matrix_s_ri_aux_kp, matrix_s, matrix_s_ri_aux, matrix_w, matrix_p_mp2, matrix_p_mp2_admm, rho, rho_xc, pw_env, ewald_env, ewald_pw, active_space, mpools, input, para_env, blacs_env, scf_control, rel_control, kinetic, qs_charges, vppl, rho_core, rho_nlcc, rho_nlcc_g, ks_env, ks_qmmm_env, wf_history, scf_env, local_particles, local_molecules, distribution_2d, dbcsr_dist, molecule_kind_set, molecule_set, subsys, cp_subsys, oce, local_rho_set, rho_atom_set, task_list, task_list_soft, rho0_atom_set, rho0_mpole, rhoz_set, rhoz_cneo_set, ecoul_1c, rho0_s_rs, rho0_s_gs, rhoz_cneo_s_rs, rhoz_cneo_s_gs, do_kpoints, has_unit_metric, requires_mo_derivs, mo_derivs, mo_loc_history, nkind, natom, nelectron_total, nelectron_spin, efield, neighbor_list_id, linres_control, xas_env, virial, cp_ddapc_env, cp_ddapc_ewald, outer_scf_history, outer_scf_ihistory, x_data, et_coupling, dftb_potential, results, se_taper, se_store_int_env, se_nddo_mpole, se_nonbond_env, admm_env, lri_env, lri_density, exstate_env, ec_env, harris_env, dispersion_env, gcp_env, vee, rho_external, external_vxc, mask, mp2_env, bs_env, kg_env, wanniercentres, atprop, ls_scf_env, do_transport, transport_env, v_hartree_rspace, s_mstruct_changed, rho_changed, potential_changed, forces_up_to_date, mscfg_env, almo_scf_env, gradient_history, variable_history, embed_pot, spin_embed_pot, polar_env, mos_last_converged, eeq, rhs, do_rixs, tb_tblite)
Get the QUICKSTEP environment.
subroutine, public prepare_gapw_den(qs_env, local_rho_set, do_rho0, kind_set_external, pw_env_sub)
...
Integrate single or product functions over a potential on a RS grid.
module that builds the second order perturbation kernel kpp1 = delta_rho|_P delta_rho|_P E drho(P1) d...
subroutine, public calc_kpp1(rho1_xc, rho1, xc_section, lrigpw, do_triplet, qs_env, p_env, calc_forces, calc_virial, virial)
...
subroutine, public kpp1_create(kpp1_env)
allocates and initializes a kpp1_env
subroutine, public kpp1_did_change(kpp1_env)
function to advise of changes either in the grids
basis types for the calculation of the perturbation of density theory.
routines that build the Kohn-Sham matrix contributions coming from local atomic densities
Definition qs_ks_atom.F:12
subroutine, public update_ks_atom(qs_env, ksmat, pmat, forces, tddft, rho_atom_external, kind_set_external, oce_external, sab_external, kscale, kintegral, kforce, fscale)
The correction to the KS matrix due to the GAPW local terms to the hartree and XC contributions is he...
Definition qs_ks_atom.F:110
basis types for the calculation of the perturbation of density theory.
subroutine, public integrate_vhg0_rspace(qs_env, v_rspace, para_env, calculate_forces, local_rho_set, local_rho_set_2nd, atener, kforce, my_pools, my_rs_descs)
...
superstucture that hold various representations of the density and keeps track of which ones are vali...
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...
routines that build the integrals of the Vxc potential calculated for the atomic density in the basis...
Definition qs_vxc_atom.F:12
subroutine, public calculate_xc_2nd_deriv_atom(rho_atom_set, rho1_atom_set, qs_env, xc_section, para_env, do_tddfpt2, do_triplet, kind_set_external)
...
represent a group ofunctional derivatives
subroutine, public xc_dset_release(derivative_set)
releases a derivative set
contains the structure
subroutine, public xc_rho_set_release(rho_set, pw_pool)
releases the given rho_set
Exchange and Correlation functional calculations.
Definition xc.F:17
subroutine, public xc_calc_2nd_deriv(v_xc, v_xc_tau, deriv_set, rho_set, rho1_r, rho1_g, tau1_r, pw_pool, xc_section, gapw, vxg, do_excitations, do_triplet, compute_virial, virial_xc)
Caller routine to calculate the second order potential in the direction of rho1_r.
Definition xc.F:1528
subroutine, public xc_prep_2nd_deriv(deriv_set, rho_set, rho_r, pw_pool, xc_section, tau_r)
Prepare objects for the calculation of the 2nd derivatives of the density functional....
Definition xc.F:5373
stores some data used in wavefunction fitting
Definition admm_types.F:120
Provides all information about an atomic kind.
type of a logger, at the moment it contains just a print level starting at which level it should be l...
stores all the informations relevant to an mpi environment
contained for different pw related things
environment for the poisson solver
Manages a pool of grids (to be used for example as tmp objects), but can also be used to instantiate ...
environment that keeps the informations and temporary val to build the kpp1 kernel matrix
Represent a qs system that is perturbed. Can calculate the linear operator and the rhs of the system ...
keeps the density in various representations, keeping track of which ones are valid.