(git:71c3ab0)
Loading...
Searching...
No Matches
qs_integrate_potential_single.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 Build up the plane wave density by collocating the primitive Gaussian
10!> functions (pgf).
11!> \par History
12!> Joost VandeVondele (02.2002)
13!> 1) rewrote collocate_pgf for increased accuracy and speed
14!> 2) collocate_core hack for PGI compiler
15!> 3) added multiple grid feature
16!> 4) new way to go over the grid
17!> Joost VandeVondele (05.2002)
18!> 1) prelim. introduction of the real space grid type
19!> JGH [30.08.02] multigrid arrays independent from potential
20!> JGH [17.07.03] distributed real space code
21!> JGH [23.11.03] refactoring and new loop ordering
22!> JGH [04.12.03] OpneMP parallelization of main loops
23!> Joost VandeVondele (12.2003)
24!> 1) modified to compute tau
25!> Joost removed incremental build feature
26!> Joost introduced map consistent
27!> Rewrote grid integration/collocation routines, [Joost VandeVondele,03.2007]
28!> \author Matthias Krack (03.04.2001)
29! **************************************************************************************************
30MODULE qs_integrate_potential_single
31 USE ao_util, ONLY: exp_radius_very_extended
32 USE atomic_kind_types, ONLY: atomic_kind_type,&
33 get_atomic_kind,&
34 get_atomic_kind_set
35 USE atprop_types, ONLY: atprop_array_init,&
36 atprop_type
37 USE basis_set_types, ONLY: get_gto_basis_set,&
38 gto_basis_set_type
39 USE cell_types, ONLY: cell_type,&
40 pbc
41 USE cp_control_types, ONLY: dft_control_type
42 USE cp_dbcsr_api, ONLY: dbcsr_get_block_p,&
43 dbcsr_type
44 USE external_potential_types, ONLY: get_potential,&
45 gth_potential_type
46 USE gaussian_gridlevels, ONLY: gaussian_gridlevel,&
47 gridlevel_info_type
48 USE grid_api, ONLY: integrate_pgf_product
49 USE kinds, ONLY: dp
50 USE lri_environment_types, ONLY: lri_kind_type
51 USE memory_utilities, ONLY: reallocate
52 USE message_passing, ONLY: mp_para_env_type
53 USE orbital_pointers, ONLY: coset,&
54 ncoset
55 USE particle_types, ONLY: particle_type
56 USE pw_env_types, ONLY: pw_env_get,&
57 pw_env_type
58 USE pw_types, ONLY: pw_r3d_rs_type
59 USE qs_environment_types, ONLY: get_qs_env,&
60 qs_environment_type
61 USE qs_force_types, ONLY: qs_force_type
62 USE qs_kind_types, ONLY: get_qs_kind,&
63 get_qs_kind_set,&
64 qs_kind_type
65 USE realspace_grid_types, ONLY: map_gaussian_here,&
66 realspace_grid_type,&
67 rs_grid_zero,&
68 transfer_pw2rs
69 USE rs_pw_interface, ONLY: potential_pw2rs
70 USE virial_types, ONLY: virial_type
71#include "./base/base_uses.f90"
72
73 IMPLICIT NONE
74
75 PRIVATE
76
77 LOGICAL, PRIVATE, PARAMETER :: debug_this_module = .false.
78
79 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'qs_integrate_potential_single'
80
81! *** Public subroutines ***
82! *** Don't include this routines directly, use the interface to
83! *** qs_integrate_potential
84
85 PUBLIC :: integrate_v_rspace_one_center, &
86 integrate_v_rspace_diagonal, &
87 integrate_v_core_rspace, &
88 integrate_v_gaussian_rspace, &
89 integrate_function, &
90 integrate_ppl_rspace, &
91 integrate_rho_nlcc
92
93CONTAINS
94
95! **************************************************************************************************
96!> \brief computes the forces/virial due to the local pseudopotential
97!> \param rho_rspace ...
98!> \param qs_env ...
99! **************************************************************************************************
100 SUBROUTINE integrate_ppl_rspace(rho_rspace, qs_env)
101 TYPE(pw_r3d_rs_type), INTENT(IN) :: rho_rspace
102 TYPE(qs_environment_type), POINTER :: qs_env
103
104 CHARACTER(len=*), PARAMETER :: routinen = 'integrate_ppl_rspace'
105
106 INTEGER :: atom_a, handle, iatom, ikind, j, lppl, &
107 n, natom_of_kind, ni, npme
108 INTEGER, DIMENSION(:), POINTER :: atom_list, cores
109 LOGICAL :: use_virial
110 REAL(kind=dp) :: alpha, eps_rho_rspace, radius
111 REAL(kind=dp), DIMENSION(3) :: force_a, force_b, ra
112 REAL(kind=dp), DIMENSION(3, 3) :: my_virial_a, my_virial_b
113 REAL(kind=dp), DIMENSION(:), POINTER :: cexp_ppl
114 REAL(kind=dp), DIMENSION(:, :), POINTER :: hab, pab
115 TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
116 TYPE(cell_type), POINTER :: cell
117 TYPE(dft_control_type), POINTER :: dft_control
118 TYPE(gth_potential_type), POINTER :: gth_potential
119 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
120 TYPE(pw_env_type), POINTER :: pw_env
121 TYPE(qs_force_type), DIMENSION(:), POINTER :: force
122 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
123 TYPE(realspace_grid_type), POINTER :: rs_v
124 TYPE(virial_type), POINTER :: virial
125
126 CALL timeset(routinen, handle)
127
128 NULLIFY (pw_env, cores)
129
130 CALL get_qs_env(qs_env=qs_env, pw_env=pw_env)
131 CALL pw_env_get(pw_env=pw_env, auxbas_rs_grid=rs_v)
132
133 CALL transfer_pw2rs(rs_v, rho_rspace)
134
135 CALL get_qs_env(qs_env=qs_env, &
136 atomic_kind_set=atomic_kind_set, &
137 qs_kind_set=qs_kind_set, &
138 cell=cell, &
139 dft_control=dft_control, &
140 particle_set=particle_set, &
141 pw_env=pw_env, &
142 force=force, virial=virial)
143
144 use_virial = virial%pv_availability .AND. (.NOT. virial%pv_numer)
145
146 eps_rho_rspace = dft_control%qs_control%eps_rho_rspace
147
148 DO ikind = 1, SIZE(atomic_kind_set)
149
150 CALL get_atomic_kind(atomic_kind_set(ikind), natom=natom_of_kind, atom_list=atom_list)
151 CALL get_qs_kind(qs_kind_set(ikind), gth_potential=gth_potential)
152
153 IF (.NOT. ASSOCIATED(gth_potential)) cycle
154 CALL get_potential(potential=gth_potential, alpha_ppl=alpha, nexp_ppl=lppl, cexp_ppl=cexp_ppl)
155
156 IF (lppl <= 0) cycle
157
158 ni = ncoset(2*lppl - 2)
159 ALLOCATE (hab(ni, 1), pab(ni, 1))
160 pab = 0._dp
161
162 CALL reallocate(cores, 1, natom_of_kind)
163 npme = 0
164 cores = 0
165
166 ! prepare core function
167 DO j = 1, lppl
168 SELECT CASE (j)
169 CASE (1)
170 pab(1, 1) = cexp_ppl(1)
171 CASE (2)
172 n = coset(2, 0, 0)
173 pab(n, 1) = cexp_ppl(2)
174 n = coset(0, 2, 0)
175 pab(n, 1) = cexp_ppl(2)
176 n = coset(0, 0, 2)
177 pab(n, 1) = cexp_ppl(2)
178 CASE (3)
179 n = coset(4, 0, 0)
180 pab(n, 1) = cexp_ppl(3)
181 n = coset(0, 4, 0)
182 pab(n, 1) = cexp_ppl(3)
183 n = coset(0, 0, 4)
184 pab(n, 1) = cexp_ppl(3)
185 n = coset(2, 2, 0)
186 pab(n, 1) = 2._dp*cexp_ppl(3)
187 n = coset(2, 0, 2)
188 pab(n, 1) = 2._dp*cexp_ppl(3)
189 n = coset(0, 2, 2)
190 pab(n, 1) = 2._dp*cexp_ppl(3)
191 CASE (4)
192 n = coset(6, 0, 0)
193 pab(n, 1) = cexp_ppl(4)
194 n = coset(0, 6, 0)
195 pab(n, 1) = cexp_ppl(4)
196 n = coset(0, 0, 6)
197 pab(n, 1) = cexp_ppl(4)
198 n = coset(4, 2, 0)
199 pab(n, 1) = 3._dp*cexp_ppl(4)
200 n = coset(4, 0, 2)
201 pab(n, 1) = 3._dp*cexp_ppl(4)
202 n = coset(2, 4, 0)
203 pab(n, 1) = 3._dp*cexp_ppl(4)
204 n = coset(2, 0, 4)
205 pab(n, 1) = 3._dp*cexp_ppl(4)
206 n = coset(0, 4, 2)
207 pab(n, 1) = 3._dp*cexp_ppl(4)
208 n = coset(0, 2, 4)
209 pab(n, 1) = 3._dp*cexp_ppl(4)
210 n = coset(2, 2, 2)
211 pab(n, 1) = 6._dp*cexp_ppl(4)
212 CASE DEFAULT
213 CALL cp_abort(__location__, &
214 "Only 1, 2, 3, 4 are supported as the "// &
215 "value of j in integrate_ppl_rspace")
216 END SELECT
217 END DO
218
219 DO iatom = 1, natom_of_kind
220 atom_a = atom_list(iatom)
221 ra(:) = pbc(particle_set(atom_a)%r, cell)
222 IF (rs_v%desc%parallel .AND. .NOT. rs_v%desc%distributed) THEN
223 ! replicated realspace grid, split the atoms up between procs
224 IF (modulo(iatom, rs_v%desc%group_size) == rs_v%desc%my_pos) THEN
225 npme = npme + 1
226 cores(npme) = iatom
227 END IF
228 ELSE
229 npme = npme + 1
230 cores(npme) = iatom
231 END IF
232 END DO
233
234 DO j = 1, npme
235
236 iatom = cores(j)
237 atom_a = atom_list(iatom)
238 ra(:) = pbc(particle_set(atom_a)%r, cell)
239 hab(:, 1) = 0.0_dp
240 force_a(:) = 0.0_dp
241 force_b(:) = 0.0_dp
242 IF (use_virial) THEN
243 my_virial_a = 0.0_dp
244 my_virial_b = 0.0_dp
245 END IF
246 ni = 2*lppl - 2
247
248 radius = exp_radius_very_extended(la_min=0, la_max=ni, lb_min=0, lb_max=0, &
249 ra=ra, rb=ra, rp=ra, &
250 zetp=alpha, eps=eps_rho_rspace, &
251 pab=pab, o1=0, o2=0, & ! without map_consistent
252 prefactor=1.0_dp, cutoff=1.0_dp)
253
254 CALL integrate_pgf_product(ni, alpha, 0, &
255 0, 0.0_dp, 0, ra, [0.0_dp, 0.0_dp, 0.0_dp], &
256 rs_v, hab, pab=pab, o1=0, o2=0, &
257 radius=radius, &
258 calculate_forces=.true., force_a=force_a, &
259 force_b=force_b, use_virial=use_virial, my_virial_a=my_virial_a, &
260 my_virial_b=my_virial_b, use_subpatch=.true., subpatch_pattern=0)
261
262 force(ikind)%gth_ppl(:, iatom) = &
263 force(ikind)%gth_ppl(:, iatom) + force_a(:)*rho_rspace%pw_grid%dvol
264
265 IF (use_virial) THEN
266 virial%pv_ppl = virial%pv_ppl + my_virial_a*rho_rspace%pw_grid%dvol
267 virial%pv_virial = virial%pv_virial + my_virial_a*rho_rspace%pw_grid%dvol
268 cpabort("Virial not debuged for CORE_PPL")
269 END IF
270 END DO
271
272 DEALLOCATE (hab, pab)
273
274 END DO
275
276 DEALLOCATE (cores)
277
278 CALL timestop(handle)
279
280 END SUBROUTINE integrate_ppl_rspace
281
282! **************************************************************************************************
283!> \brief computes the forces/virial due to the nlcc pseudopotential
284!> \param rho_rspace ...
285!> \param qs_env ...
286! **************************************************************************************************
287 SUBROUTINE integrate_rho_nlcc(rho_rspace, qs_env)
288 TYPE(pw_r3d_rs_type), INTENT(IN) :: rho_rspace
289 TYPE(qs_environment_type), POINTER :: qs_env
290
291 CHARACTER(len=*), PARAMETER :: routinen = 'integrate_rho_nlcc'
292
293 INTEGER :: atom_a, handle, iatom, iexp_nlcc, ikind, &
294 ithread, j, n, natom, nc, nexp_nlcc, &
295 ni, npme, nthread
296 INTEGER, DIMENSION(:), POINTER :: atom_list, cores, nct_nlcc
297 LOGICAL :: nlcc, use_virial
298 REAL(kind=dp) :: alpha, eps_rho_rspace, radius
299 REAL(kind=dp), DIMENSION(3) :: force_a, force_b, ra
300 REAL(kind=dp), DIMENSION(3, 3) :: my_virial_a, my_virial_b
301 REAL(kind=dp), DIMENSION(:), POINTER :: alpha_nlcc
302 REAL(kind=dp), DIMENSION(:, :), POINTER :: cval_nlcc, hab, pab
303 TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
304 TYPE(cell_type), POINTER :: cell
305 TYPE(dft_control_type), POINTER :: dft_control
306 TYPE(gth_potential_type), POINTER :: gth_potential
307 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
308 TYPE(pw_env_type), POINTER :: pw_env
309 TYPE(qs_force_type), DIMENSION(:), POINTER :: force
310 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
311 TYPE(realspace_grid_type), POINTER :: rs_v
312 TYPE(virial_type), POINTER :: virial
313
314 CALL timeset(routinen, handle)
315
316 NULLIFY (pw_env, cores)
317
318 CALL get_qs_env(qs_env=qs_env, pw_env=pw_env)
319 CALL pw_env_get(pw_env=pw_env, auxbas_rs_grid=rs_v)
320
321 CALL transfer_pw2rs(rs_v, rho_rspace)
322
323 CALL get_qs_env(qs_env=qs_env, &
324 atomic_kind_set=atomic_kind_set, &
325 qs_kind_set=qs_kind_set, &
326 cell=cell, &
327 dft_control=dft_control, &
328 particle_set=particle_set, &
329 pw_env=pw_env, &
330 force=force, virial=virial)
331
332 use_virial = virial%pv_availability .AND. (.NOT. virial%pv_numer)
333
334 eps_rho_rspace = dft_control%qs_control%eps_rho_rspace
335
336 DO ikind = 1, SIZE(atomic_kind_set)
337
338 CALL get_atomic_kind(atomic_kind_set(ikind), natom=natom, atom_list=atom_list)
339 CALL get_qs_kind(qs_kind_set(ikind), gth_potential=gth_potential)
340
341 IF (.NOT. ASSOCIATED(gth_potential)) cycle
342 CALL get_potential(potential=gth_potential, nlcc_present=nlcc, nexp_nlcc=nexp_nlcc, &
343 alpha_nlcc=alpha_nlcc, nct_nlcc=nct_nlcc, cval_nlcc=cval_nlcc)
344
345 IF (.NOT. nlcc) cycle
346
347 DO iexp_nlcc = 1, nexp_nlcc
348
349 alpha = alpha_nlcc(iexp_nlcc)
350 nc = nct_nlcc(iexp_nlcc)
351
352 ni = ncoset(2*nc - 2)
353
354 nthread = 1
355 ithread = 0
356
357 ALLOCATE (hab(ni, 1), pab(ni, 1))
358 pab = 0._dp
359
360 CALL reallocate(cores, 1, natom)
361 npme = 0
362 cores = 0
363
364 ! prepare core function
365 DO j = 1, nc
366 SELECT CASE (j)
367 CASE (1)
368 pab(1, 1) = cval_nlcc(1, iexp_nlcc)
369 CASE (2)
370 n = coset(2, 0, 0)
371 pab(n, 1) = cval_nlcc(2, iexp_nlcc)/alpha**2
372 n = coset(0, 2, 0)
373 pab(n, 1) = cval_nlcc(2, iexp_nlcc)/alpha**2
374 n = coset(0, 0, 2)
375 pab(n, 1) = cval_nlcc(2, iexp_nlcc)/alpha**2
376 CASE (3)
377 n = coset(4, 0, 0)
378 pab(n, 1) = cval_nlcc(3, iexp_nlcc)/alpha**4
379 n = coset(0, 4, 0)
380 pab(n, 1) = cval_nlcc(3, iexp_nlcc)/alpha**4
381 n = coset(0, 0, 4)
382 pab(n, 1) = cval_nlcc(3, iexp_nlcc)/alpha**4
383 n = coset(2, 2, 0)
384 pab(n, 1) = 2._dp*cval_nlcc(3, iexp_nlcc)/alpha**4
385 n = coset(2, 0, 2)
386 pab(n, 1) = 2._dp*cval_nlcc(3, iexp_nlcc)/alpha**4
387 n = coset(0, 2, 2)
388 pab(n, 1) = 2._dp*cval_nlcc(3, iexp_nlcc)/alpha**4
389 CASE (4)
390 n = coset(6, 0, 0)
391 pab(n, 1) = cval_nlcc(4, iexp_nlcc)/alpha**6
392 n = coset(0, 6, 0)
393 pab(n, 1) = cval_nlcc(4, iexp_nlcc)/alpha**6
394 n = coset(0, 0, 6)
395 pab(n, 1) = cval_nlcc(4, iexp_nlcc)/alpha**6
396 n = coset(4, 2, 0)
397 pab(n, 1) = 3._dp*cval_nlcc(4, iexp_nlcc)/alpha**6
398 n = coset(4, 0, 2)
399 pab(n, 1) = 3._dp*cval_nlcc(4, iexp_nlcc)/alpha**6
400 n = coset(2, 4, 0)
401 pab(n, 1) = 3._dp*cval_nlcc(4, iexp_nlcc)/alpha**6
402 n = coset(2, 0, 4)
403 pab(n, 1) = 3._dp*cval_nlcc(4, iexp_nlcc)/alpha**6
404 n = coset(0, 4, 2)
405 pab(n, 1) = 3._dp*cval_nlcc(4, iexp_nlcc)/alpha**6
406 n = coset(0, 2, 4)
407 pab(n, 1) = 3._dp*cval_nlcc(4, iexp_nlcc)/alpha**6
408 n = coset(2, 2, 2)
409 pab(n, 1) = 6._dp*cval_nlcc(4, iexp_nlcc)/alpha**6
410 CASE DEFAULT
411 CALL cp_abort(__location__, &
412 "Only 1, 2, 3, 4 are supported as the "// &
413 "value of j in integrate_rho_nlcc")
414 END SELECT
415 END DO
416 IF (dft_control%nspins == 2) pab = pab*0.5_dp
417
418 DO iatom = 1, natom
419 atom_a = atom_list(iatom)
420 ra(:) = pbc(particle_set(atom_a)%r, cell)
421 IF (rs_v%desc%parallel .AND. .NOT. rs_v%desc%distributed) THEN
422 ! replicated realspace grid, split the atoms up between procs
423 IF (modulo(iatom, rs_v%desc%group_size) == rs_v%desc%my_pos) THEN
424 npme = npme + 1
425 cores(npme) = iatom
426 END IF
427 ELSE
428 npme = npme + 1
429 cores(npme) = iatom
430 END IF
431 END DO
432
433 DO j = 1, npme
434
435 iatom = cores(j)
436 atom_a = atom_list(iatom)
437 ra(:) = pbc(particle_set(atom_a)%r, cell)
438 hab(:, 1) = 0.0_dp
439 force_a(:) = 0.0_dp
440 force_b(:) = 0.0_dp
441 IF (use_virial) THEN
442 my_virial_a = 0.0_dp
443 my_virial_b = 0.0_dp
444 END IF
445 ni = 2*nc - 2
446
447 radius = exp_radius_very_extended(la_min=0, la_max=ni, lb_min=0, lb_max=0, &
448 ra=ra, rb=ra, rp=ra, &
449 zetp=1/(2*alpha**2), eps=eps_rho_rspace, &
450 pab=pab, o1=0, o2=0, & ! without map_consistent
451 prefactor=1.0_dp, cutoff=1.0_dp)
452
453 CALL integrate_pgf_product(ni, 1/(2*alpha**2), 0, &
454 0, 0.0_dp, 0, ra, [0.0_dp, 0.0_dp, 0.0_dp], &
455 rs_v, hab, pab=pab, o1=0, o2=0, &
456 radius=radius, &
457 calculate_forces=.true., force_a=force_a, &
458 force_b=force_b, use_virial=use_virial, my_virial_a=my_virial_a, &
459 my_virial_b=my_virial_b, use_subpatch=.true., subpatch_pattern=0)
460
461 force(ikind)%gth_nlcc(:, iatom) = &
462 force(ikind)%gth_nlcc(:, iatom) + force_a(:)*rho_rspace%pw_grid%dvol
463
464 IF (use_virial) THEN
465 virial%pv_nlcc = virial%pv_nlcc + my_virial_a*rho_rspace%pw_grid%dvol
466 virial%pv_virial = virial%pv_virial + my_virial_a*rho_rspace%pw_grid%dvol
467 END IF
468 END DO
469
470 DEALLOCATE (hab, pab)
471
472 END DO
473
474 END DO
475
476 DEALLOCATE (cores)
477
478 CALL timestop(handle)
479
480 END SUBROUTINE integrate_rho_nlcc
481
482! **************************************************************************************************
483!> \brief computes the forces/virial due to the ionic cores with a potential on
484!> grid
485!> \param v_rspace ...
486!> \param qs_env ...
487!> \param atecc ...
488! **************************************************************************************************
489 SUBROUTINE integrate_v_core_rspace(v_rspace, qs_env, atecc)
490 TYPE(pw_r3d_rs_type), INTENT(IN) :: v_rspace
491 TYPE(qs_environment_type), POINTER :: qs_env
492 REAL(kind=dp), DIMENSION(:), OPTIONAL :: atecc
493
494 CHARACTER(len=*), PARAMETER :: routinen = 'integrate_v_core_rspace'
495
496 INTEGER :: atom_a, handle, iatom, ikind, j, natom, &
497 natom_of_kind, npme
498 INTEGER, DIMENSION(:), POINTER :: atom_list, cores
499 LOGICAL :: paw_atom, skip_fcore, use_virial
500 REAL(kind=dp) :: alpha_core_charge, ccore_charge, &
501 eps_rho_rspace, radius
502 REAL(kind=dp), DIMENSION(3) :: force_a, force_b, ra
503 REAL(kind=dp), DIMENSION(3, 3) :: my_virial_a, my_virial_b
504 REAL(kind=dp), DIMENSION(:, :), POINTER :: hab, pab
505 TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
506 TYPE(atprop_type), POINTER :: atprop
507 TYPE(cell_type), POINTER :: cell
508 TYPE(dft_control_type), POINTER :: dft_control
509 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
510 TYPE(pw_env_type), POINTER :: pw_env
511 TYPE(qs_force_type), DIMENSION(:), POINTER :: force
512 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
513 TYPE(realspace_grid_type), POINTER :: rs_v
514 TYPE(virial_type), POINTER :: virial
515
516 CALL timeset(routinen, handle)
517 NULLIFY (virial, force, atprop, dft_control)
518
519 CALL get_qs_env(qs_env=qs_env, dft_control=dft_control)
520
521 !If gapw, check for gpw kinds
522 skip_fcore = .false.
523 IF (dft_control%qs_control%gapw) THEN
524 IF (.NOT. dft_control%qs_control%gapw_control%nopaw_as_gpw) skip_fcore = .true.
525 END IF
526
527 IF (.NOT. skip_fcore) THEN
528 NULLIFY (pw_env)
529 ALLOCATE (cores(1))
530 ALLOCATE (hab(1, 1))
531 ALLOCATE (pab(1, 1))
532
533 CALL get_qs_env(qs_env=qs_env, pw_env=pw_env)
534 CALL pw_env_get(pw_env=pw_env, auxbas_rs_grid=rs_v)
535
536 CALL transfer_pw2rs(rs_v, v_rspace)
537
538 CALL get_qs_env(qs_env=qs_env, &
539 atomic_kind_set=atomic_kind_set, &
540 qs_kind_set=qs_kind_set, &
541 cell=cell, &
542 dft_control=dft_control, &
543 particle_set=particle_set, &
544 pw_env=pw_env, &
545 force=force, &
546 virial=virial, &
547 atprop=atprop)
548
549 ! atomic energy contributions
550 natom = SIZE(particle_set)
551 IF (ASSOCIATED(atprop)) THEN
552 CALL atprop_array_init(atprop%ateb, natom)
553 END IF
554
555 use_virial = virial%pv_availability .AND. (.NOT. virial%pv_numer)
556
557 eps_rho_rspace = dft_control%qs_control%eps_rho_rspace
558
559 DO ikind = 1, SIZE(atomic_kind_set)
560
561 CALL get_atomic_kind(atomic_kind_set(ikind), natom=natom_of_kind, atom_list=atom_list)
562 CALL get_qs_kind(qs_kind_set(ikind), paw_atom=paw_atom)
563
564 IF (dft_control%qs_control%gapw .AND. paw_atom) cycle
565
566 CALL get_qs_kind(qs_kind_set(ikind), alpha_core_charge=alpha_core_charge, &
567 ccore_charge=ccore_charge)
568
569 pab(1, 1) = -ccore_charge
570 IF (alpha_core_charge == 0.0_dp .OR. pab(1, 1) == 0.0_dp) cycle
571
572 CALL reallocate(cores, 1, natom_of_kind)
573 npme = 0
574 cores = 0
575
576 DO iatom = 1, natom_of_kind
577 atom_a = atom_list(iatom)
578 ra(:) = pbc(particle_set(atom_a)%r, cell)
579 IF (rs_v%desc%parallel .AND. .NOT. rs_v%desc%distributed) THEN
580 ! replicated realspace grid, split the atoms up between procs
581 IF (modulo(iatom, rs_v%desc%group_size) == rs_v%desc%my_pos) THEN
582 npme = npme + 1
583 cores(npme) = iatom
584 END IF
585 ELSE
586 npme = npme + 1
587 cores(npme) = iatom
588 END IF
589 END DO
590
591 DO j = 1, npme
592
593 iatom = cores(j)
594 atom_a = atom_list(iatom)
595 ra(:) = pbc(particle_set(atom_a)%r, cell)
596 hab(1, 1) = 0.0_dp
597 force_a(:) = 0.0_dp
598 force_b(:) = 0.0_dp
599 IF (use_virial) THEN
600 my_virial_a = 0.0_dp
601 my_virial_b = 0.0_dp
602 END IF
603
604 radius = exp_radius_very_extended(la_min=0, la_max=0, lb_min=0, lb_max=0, &
605 ra=ra, rb=ra, rp=ra, &
606 zetp=alpha_core_charge, eps=eps_rho_rspace, &
607 pab=pab, o1=0, o2=0, & ! without map_consistent
608 prefactor=1.0_dp, cutoff=1.0_dp)
609
610 CALL integrate_pgf_product(0, alpha_core_charge, 0, &
611 0, 0.0_dp, 0, ra, [0.0_dp, 0.0_dp, 0.0_dp], &
612 rs_v, hab, pab=pab, o1=0, o2=0, &
613 radius=radius, &
614 calculate_forces=.true., force_a=force_a, &
615 force_b=force_b, use_virial=use_virial, my_virial_a=my_virial_a, &
616 my_virial_b=my_virial_b, use_subpatch=.true., subpatch_pattern=0)
617
618 IF (ASSOCIATED(force)) THEN
619 force(ikind)%rho_core(:, iatom) = force(ikind)%rho_core(:, iatom) + force_a(:)
620 END IF
621 IF (use_virial) THEN
622 virial%pv_ehartree = virial%pv_ehartree + my_virial_a
623 virial%pv_virial = virial%pv_virial + my_virial_a
624 END IF
625 IF (ASSOCIATED(atprop)) THEN
626 atprop%ateb(atom_a) = atprop%ateb(atom_a) + 0.5_dp*hab(1, 1)*pab(1, 1)
627 END IF
628 IF (PRESENT(atecc)) THEN
629 atecc(atom_a) = atecc(atom_a) + 0.5_dp*hab(1, 1)*pab(1, 1)
630 END IF
631
632 END DO
633
634 END DO
635
636 DEALLOCATE (hab, pab, cores)
637
638 END IF
639
640 CALL timestop(handle)
641
642 END SUBROUTINE integrate_v_core_rspace
643
644! **************************************************************************************************
645!> \brief computes the overlap of a set of Gaussians with a potential on grid
646!> \param v_rspace ...
647!> \param qs_env ...
648!> \param alpha ...
649!> \param ccore ...
650!> \param atecc ...
651! **************************************************************************************************
652 SUBROUTINE integrate_v_gaussian_rspace(v_rspace, qs_env, alpha, ccore, atecc)
653 TYPE(pw_r3d_rs_type), INTENT(IN) :: v_rspace
654 TYPE(qs_environment_type), POINTER :: qs_env
655 REAL(kind=dp), DIMENSION(:), INTENT(IN) :: alpha, ccore
656 REAL(kind=dp), DIMENSION(:) :: atecc
657
658 CHARACTER(len=*), PARAMETER :: routinen = 'integrate_v_gaussian_rspace'
659
660 INTEGER :: atom_a, handle, iatom, ikind, j, natom, &
661 natom_of_kind, npme
662 INTEGER, DIMENSION(:), POINTER :: atom_list, cores
663 REAL(kind=dp) :: alpha_core_charge, eps_rho_rspace, radius
664 REAL(kind=dp), DIMENSION(3) :: ra
665 REAL(kind=dp), DIMENSION(:, :), POINTER :: hab, pab
666 TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
667 TYPE(cell_type), POINTER :: cell
668 TYPE(dft_control_type), POINTER :: dft_control
669 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
670 TYPE(pw_env_type), POINTER :: pw_env
671 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
672 TYPE(realspace_grid_type), POINTER :: rs_v
673
674 CALL timeset(routinen, handle)
675
676 CALL get_qs_env(qs_env=qs_env, dft_control=dft_control)
677
678 !If gapw, check for gpw kinds
679 cpassert(.NOT. dft_control%qs_control%gapw)
680
681 NULLIFY (pw_env)
682 ALLOCATE (cores(1))
683 ALLOCATE (hab(1, 1))
684 ALLOCATE (pab(1, 1))
685
686 CALL get_qs_env(qs_env=qs_env, pw_env=pw_env)
687 CALL pw_env_get(pw_env=pw_env, auxbas_rs_grid=rs_v)
688
689 CALL transfer_pw2rs(rs_v, v_rspace)
690
691 CALL get_qs_env(qs_env=qs_env, &
692 atomic_kind_set=atomic_kind_set, &
693 qs_kind_set=qs_kind_set, &
694 cell=cell, &
695 dft_control=dft_control, &
696 particle_set=particle_set, &
697 pw_env=pw_env)
698
699 ! atomic energy contributions
700 natom = SIZE(particle_set)
701 eps_rho_rspace = dft_control%qs_control%eps_rho_rspace
702
703 DO ikind = 1, SIZE(atomic_kind_set)
704
705 CALL get_atomic_kind(atomic_kind_set(ikind), natom=natom_of_kind, atom_list=atom_list)
706 pab(1, 1) = -ccore(ikind)
707 alpha_core_charge = alpha(ikind)
708 IF (alpha_core_charge == 0.0_dp .OR. pab(1, 1) == 0.0_dp) cycle
709
710 CALL reallocate(cores, 1, natom_of_kind)
711 npme = 0
712 cores = 0
713
714 DO iatom = 1, natom_of_kind
715 atom_a = atom_list(iatom)
716 ra(:) = pbc(particle_set(atom_a)%r, cell)
717 IF (rs_v%desc%parallel .AND. .NOT. rs_v%desc%distributed) THEN
718 ! replicated realspace grid, split the atoms up between procs
719 IF (modulo(iatom, rs_v%desc%group_size) == rs_v%desc%my_pos) THEN
720 npme = npme + 1
721 cores(npme) = iatom
722 END IF
723 ELSE
724 npme = npme + 1
725 cores(npme) = iatom
726 END IF
727 END DO
728
729 DO j = 1, npme
730
731 iatom = cores(j)
732 atom_a = atom_list(iatom)
733 ra(:) = pbc(particle_set(atom_a)%r, cell)
734 hab(1, 1) = 0.0_dp
735
736 radius = exp_radius_very_extended(la_min=0, la_max=0, lb_min=0, lb_max=0, &
737 ra=ra, rb=ra, rp=ra, &
738 zetp=alpha_core_charge, eps=eps_rho_rspace, &
739 pab=pab, o1=0, o2=0, & ! without map_consistent
740 prefactor=1.0_dp, cutoff=1.0_dp)
741
742 CALL integrate_pgf_product(0, alpha_core_charge, 0, &
743 0, 0.0_dp, 0, ra, [0.0_dp, 0.0_dp, 0.0_dp], &
744 rs_v, hab, pab=pab, o1=0, o2=0, &
745 radius=radius, calculate_forces=.false., &
746 use_subpatch=.true., subpatch_pattern=0)
747 atecc(atom_a) = atecc(atom_a) + 0.5_dp*hab(1, 1)*pab(1, 1)
748
749 END DO
750
751 END DO
752
753 DEALLOCATE (hab, pab, cores)
754
755 CALL timestop(handle)
756
757 END SUBROUTINE integrate_v_gaussian_rspace
758! **************************************************************************************************
759!> \brief computes integrals of product of v_rspace times a one-center function
760!> required for LRIGPW
761!> \param v_rspace ...
762!> \param qs_env ...
763!> \param int_res ...
764!> \param calculate_forces ...
765!> \param basis_type ...
766!> \param atomlist ...
767!> \author Dorothea Golze
768! **************************************************************************************************
769 SUBROUTINE integrate_v_rspace_one_center(v_rspace, qs_env, int_res, &
770 calculate_forces, basis_type, atomlist)
771 TYPE(pw_r3d_rs_type), INTENT(IN) :: v_rspace
772 TYPE(qs_environment_type), POINTER :: qs_env
773 TYPE(lri_kind_type), DIMENSION(:), POINTER :: int_res
774 LOGICAL, INTENT(IN) :: calculate_forces
775 CHARACTER(len=*), INTENT(IN) :: basis_type
776 INTEGER, DIMENSION(:), OPTIONAL :: atomlist
777
778 CHARACTER(len=*), PARAMETER :: routinen = 'integrate_v_rspace_one_center'
779
780 INTEGER :: atom_a, group_size, handle, i, iatom, igrid_level, ikind, ipgf, iset, m1, &
781 max_npgf, maxco, maxsgf_set, my_pos, na1, natom_of_kind, ncoa, nkind, nseta, offset, sgfa
782 INTEGER, DIMENSION(:), POINTER :: atom_list, la_max, la_min, npgfa, &
783 nsgf_seta
784 INTEGER, DIMENSION(:, :), POINTER :: first_sgfa
785 LOGICAL :: use_virial
786 LOGICAL, ALLOCATABLE, DIMENSION(:) :: map_it
787 REAL(kind=dp) :: eps_rho_rspace, radius
788 REAL(kind=dp), DIMENSION(3) :: force_a, force_b, ra
789 REAL(kind=dp), DIMENSION(3, 3) :: my_virial_a, my_virial_b
790 REAL(kind=dp), DIMENSION(:), POINTER :: set_radius_a
791 REAL(kind=dp), DIMENSION(:, :), POINTER :: hab, pab, rpgfa, sphi_a, work_f, work_i, &
792 zeta
793 TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
794 TYPE(cell_type), POINTER :: cell
795 TYPE(dft_control_type), POINTER :: dft_control
796 TYPE(gridlevel_info_type), POINTER :: gridlevel_info
797 TYPE(gto_basis_set_type), POINTER :: lri_basis_set
798 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
799 TYPE(pw_env_type), POINTER :: pw_env
800 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
801 TYPE(realspace_grid_type), DIMENSION(:), POINTER :: rs_v
802 TYPE(realspace_grid_type), POINTER :: rs_grid
803 TYPE(virial_type), POINTER :: virial
804
805 CALL timeset(routinen, handle)
806
807 NULLIFY (atomic_kind_set, qs_kind_set, atom_list, cell, dft_control, &
808 first_sgfa, gridlevel_info, hab, la_max, la_min, lri_basis_set, &
809 npgfa, nsgf_seta, pab, particle_set, pw_env, rpgfa, &
810 rs_grid, rs_v, virial, set_radius_a, sphi_a, work_f, &
811 work_i, zeta)
812
813 CALL get_qs_env(qs_env=qs_env, pw_env=pw_env)
814
815 CALL pw_env_get(pw_env, rs_grids=rs_v)
816 DO i = 1, SIZE(rs_v)
817 CALL rs_grid_zero(rs_v(i))
818 END DO
819
820 gridlevel_info => pw_env%gridlevel_info
821
822 CALL potential_pw2rs(rs_v, v_rspace, pw_env)
823
824 CALL get_qs_env(qs_env=qs_env, &
825 atomic_kind_set=atomic_kind_set, &
826 qs_kind_set=qs_kind_set, &
827 cell=cell, &
828 dft_control=dft_control, &
829 nkind=nkind, &
830 particle_set=particle_set, &
831 pw_env=pw_env, &
832 virial=virial)
833
834 use_virial = virial%pv_availability .AND. (.NOT. virial%pv_numer)
835
836 eps_rho_rspace = dft_control%qs_control%eps_rho_rspace
837
838 offset = 0
839 my_pos = v_rspace%pw_grid%para%group%mepos
840 group_size = v_rspace%pw_grid%para%group%num_pe
841
842 DO ikind = 1, nkind
843
844 CALL get_atomic_kind(atomic_kind_set(ikind), natom=natom_of_kind, atom_list=atom_list)
845 CALL get_qs_kind(qs_kind_set(ikind), basis_set=lri_basis_set, basis_type=basis_type)
846 CALL get_gto_basis_set(gto_basis_set=lri_basis_set, &
847 first_sgf=first_sgfa, &
848 lmax=la_max, &
849 lmin=la_min, &
850 maxco=maxco, &
851 maxsgf_set=maxsgf_set, &
852 npgf=npgfa, &
853 nset=nseta, &
854 nsgf_set=nsgf_seta, &
855 pgf_radius=rpgfa, &
856 set_radius=set_radius_a, &
857 sphi=sphi_a, &
858 zet=zeta)
859
860 ALLOCATE (hab(maxco, 1), pab(maxco, 1))
861 hab = 0._dp
862 pab(:, 1) = 0._dp
863 max_npgf = maxval(npgfa(1:nseta))
864 ALLOCATE (map_it(max_npgf))
865 ALLOCATE (work_i(maxsgf_set, 1))
866 IF (calculate_forces) ALLOCATE (work_f(maxsgf_set, 1))
867
868 DO iatom = 1, natom_of_kind
869
870 atom_a = atom_list(iatom)
871 IF (PRESENT(atomlist)) THEN
872 IF (atomlist(atom_a) == 0) cycle
873 END IF
874 ra(:) = pbc(particle_set(atom_a)%r, cell)
875 force_a(:) = 0._dp
876 force_b(:) = 0._dp
877 my_virial_a(:, :) = 0._dp
878 my_virial_b(:, :) = 0._dp
879
880 DO iset = 1, nseta
881 !
882 map_it = .false.
883 DO ipgf = 1, npgfa(iset)
884 igrid_level = gaussian_gridlevel(gridlevel_info, zeta(ipgf, iset))
885 rs_grid => rs_v(igrid_level)
886 map_it(ipgf) = map_gaussian_here(rs_grid, cell%h_inv, ra, offset, group_size, my_pos)
887 END DO
888 offset = offset + 1
889 !
890 IF (any(map_it(1:npgfa(iset)))) THEN
891 sgfa = first_sgfa(1, iset)
892 ncoa = npgfa(iset)*ncoset(la_max(iset))
893 hab(:, 1) = 0._dp
894 work_i(1:nsgf_seta(iset), 1) = 0.0_dp
895
896 ! get fit coefficients for forces
897 IF (calculate_forces) THEN
898 m1 = sgfa + nsgf_seta(iset) - 1
899 work_f(1:nsgf_seta(iset), 1) = int_res(ikind)%acoef(iatom, sgfa:m1)
900 CALL dgemm("N", "N", ncoa, 1, nsgf_seta(iset), 1.0_dp, sphi_a(1, sgfa), &
901 SIZE(sphi_a, 1), work_f(1, 1), SIZE(work_f, 1), 0.0_dp, pab(1, 1), &
902 SIZE(pab, 1))
903 END IF
904
905 DO ipgf = 1, npgfa(iset)
906 na1 = (ipgf - 1)*ncoset(la_max(iset))
907 igrid_level = gaussian_gridlevel(gridlevel_info, zeta(ipgf, iset))
908 rs_grid => rs_v(igrid_level)
909
910 radius = exp_radius_very_extended(la_min=la_min(iset), la_max=la_max(iset), &
911 lb_min=0, lb_max=0, ra=ra, rb=ra, rp=ra, &
912 zetp=zeta(ipgf, iset), eps=eps_rho_rspace, &
913 prefactor=1.0_dp, cutoff=1.0_dp)
914
915 IF (map_it(ipgf)) THEN
916 IF (.NOT. calculate_forces) THEN
917 CALL integrate_pgf_product(la_max=la_max(iset), &
918 zeta=zeta(ipgf, iset), la_min=la_min(iset), &
919 lb_max=0, zetb=0.0_dp, lb_min=0, &
920 ra=ra, rab=[0.0_dp, 0.0_dp, 0.0_dp], &
921 rsgrid=rs_grid, &
922 hab=hab, o1=na1, o2=0, radius=radius, &
923 calculate_forces=calculate_forces)
924 ELSE
925 CALL integrate_pgf_product(la_max=la_max(iset), &
926 zeta=zeta(ipgf, iset), la_min=la_min(iset), &
927 lb_max=0, zetb=0.0_dp, lb_min=0, &
928 ra=ra, rab=[0.0_dp, 0.0_dp, 0.0_dp], &
929 rsgrid=rs_grid, &
930 hab=hab, pab=pab, o1=na1, o2=0, radius=radius, &
931 calculate_forces=calculate_forces, &
932 force_a=force_a, force_b=force_b, &
933 use_virial=use_virial, &
934 my_virial_a=my_virial_a, my_virial_b=my_virial_b)
935 END IF
936 END IF
937 END DO
938 ! contract hab
939 CALL dgemm("T", "N", nsgf_seta(iset), 1, ncoa, 1.0_dp, sphi_a(1, sgfa), &
940 SIZE(sphi_a, 1), hab(1, 1), SIZE(hab, 1), 0.0_dp, work_i(1, 1), SIZE(work_i, 1))
941
942 int_res(ikind)%v_int(iatom, sgfa:sgfa - 1 + nsgf_seta(iset)) = &
943 int_res(ikind)%v_int(iatom, sgfa:sgfa - 1 + nsgf_seta(iset)) + work_i(1:nsgf_seta(iset), 1)
944 END IF
945 END DO
946 !
947 IF (calculate_forces) THEN
948 int_res(ikind)%v_dfdr(iatom, :) = int_res(ikind)%v_dfdr(iatom, :) + force_a(:)
949 IF (use_virial) THEN
950 virial%pv_lrigpw = virial%pv_lrigpw + my_virial_a
951 virial%pv_virial = virial%pv_virial + my_virial_a
952 END IF
953 END IF
954
955 END DO
956
957 IF (calculate_forces) DEALLOCATE (work_f)
958 DEALLOCATE (work_i, map_it)
959 DEALLOCATE (hab, pab)
960 END DO
961
962 CALL timestop(handle)
963
964 END SUBROUTINE integrate_v_rspace_one_center
965
966! **************************************************************************************************
967!> \brief computes integrals of product of v_rspace times the diagonal block basis functions
968!> required for LRIGPW with exact 1c terms
969!> \param v_rspace ...
970!> \param ksmat ...
971!> \param pmat ...
972!> \param qs_env ...
973!> \param calculate_forces ...
974!> \param basis_type ...
975!> \author JGH
976! **************************************************************************************************
977 SUBROUTINE integrate_v_rspace_diagonal(v_rspace, ksmat, pmat, qs_env, calculate_forces, basis_type)
978 TYPE(pw_r3d_rs_type), INTENT(IN) :: v_rspace
979 TYPE(dbcsr_type), INTENT(INOUT) :: ksmat, pmat
980 TYPE(qs_environment_type), POINTER :: qs_env
981 LOGICAL, INTENT(IN) :: calculate_forces
982 CHARACTER(len=*), INTENT(IN) :: basis_type
983
984 CHARACTER(len=*), PARAMETER :: routinen = 'integrate_v_rspace_diagonal'
985
986 INTEGER :: atom_a, group_size, handle, iatom, igrid_level, ikind, ipgf, iset, jpgf, jset, &
987 m1, maxco, maxsgf_set, my_pos, na1, na2, natom_of_kind, nb1, nb2, ncoa, ncob, nkind, &
988 nseta, nsgfa, offset, sgfa, sgfb
989 INTEGER, DIMENSION(:), POINTER :: atom_list, la_max, la_min, npgfa, &
990 nsgf_seta
991 INTEGER, DIMENSION(:, :), POINTER :: first_sgfa
992 LOGICAL :: found, use_virial
993 LOGICAL, ALLOCATABLE, DIMENSION(:, :) :: map_it2
994 REAL(kind=dp) :: eps_rho_rspace, radius, zetp
995 REAL(kind=dp), DIMENSION(3) :: force_a, force_b, ra
996 REAL(kind=dp), DIMENSION(3, 3) :: my_virial_a, my_virial_b
997 REAL(kind=dp), DIMENSION(:), POINTER :: set_radius_a
998 REAL(kind=dp), DIMENSION(:, :), POINTER :: h_block, hab, hmat, p_block, pab, pblk, &
999 rpgfa, sphi_a, work, zeta
1000 TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
1001 TYPE(cell_type), POINTER :: cell
1002 TYPE(dft_control_type), POINTER :: dft_control
1003 TYPE(gridlevel_info_type), POINTER :: gridlevel_info
1004 TYPE(gto_basis_set_type), POINTER :: orb_basis_set
1005 TYPE(mp_para_env_type), POINTER :: para_env
1006 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
1007 TYPE(pw_env_type), POINTER :: pw_env
1008 TYPE(qs_force_type), DIMENSION(:), POINTER :: force
1009 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
1010 TYPE(realspace_grid_type), DIMENSION(:), POINTER :: rs_v
1011 TYPE(realspace_grid_type), POINTER :: rs_grid
1012 TYPE(virial_type), POINTER :: virial
1013
1014 CALL timeset(routinen, handle)
1015
1016 CALL get_qs_env(qs_env=qs_env, pw_env=pw_env)
1017 CALL pw_env_get(pw_env, rs_grids=rs_v)
1018 CALL potential_pw2rs(rs_v, v_rspace, pw_env)
1019
1020 gridlevel_info => pw_env%gridlevel_info
1021
1022 CALL get_qs_env(qs_env=qs_env, &
1023 atomic_kind_set=atomic_kind_set, &
1024 qs_kind_set=qs_kind_set, &
1025 cell=cell, &
1026 dft_control=dft_control, &
1027 nkind=nkind, &
1028 particle_set=particle_set, &
1029 force=force, &
1030 virial=virial, &
1031 para_env=para_env)
1032
1033 eps_rho_rspace = dft_control%qs_control%eps_rho_rspace
1034 use_virial = virial%pv_availability .AND. (.NOT. virial%pv_numer)
1035
1036 offset = 0
1037 my_pos = v_rspace%pw_grid%para%group%mepos
1038 group_size = v_rspace%pw_grid%para%group%num_pe
1039
1040 DO ikind = 1, nkind
1041
1042 CALL get_atomic_kind(atomic_kind_set(ikind), natom=natom_of_kind, atom_list=atom_list)
1043 CALL get_qs_kind(qs_kind_set(ikind), basis_set=orb_basis_set, basis_type=basis_type)
1044 CALL get_gto_basis_set(gto_basis_set=orb_basis_set, &
1045 lmax=la_max, lmin=la_min, maxco=maxco, maxsgf_set=maxsgf_set, &
1046 npgf=npgfa, nset=nseta, nsgf_set=nsgf_seta, nsgf=nsgfa, &
1047 first_sgf=first_sgfa, pgf_radius=rpgfa, set_radius=set_radius_a, &
1048 sphi=sphi_a, zet=zeta)
1049
1050 ALLOCATE (hab(maxco, maxco), work(maxco, maxsgf_set), hmat(nsgfa, nsgfa))
1051 IF (calculate_forces) ALLOCATE (pab(maxco, maxco), pblk(nsgfa, nsgfa))
1052
1053 DO iatom = 1, natom_of_kind
1054 atom_a = atom_list(iatom)
1055 ra(:) = pbc(particle_set(atom_a)%r, cell)
1056 hmat = 0.0_dp
1057 IF (calculate_forces) THEN
1058 CALL dbcsr_get_block_p(matrix=pmat, row=atom_a, col=atom_a, block=p_block, found=found)
1059 IF (found) THEN
1060 pblk(1:nsgfa, 1:nsgfa) = p_block(1:nsgfa, 1:nsgfa)
1061 ELSE
1062 pblk = 0.0_dp
1063 END IF
1064 CALL para_env%sum(pblk)
1065 force_a(:) = 0._dp
1066 force_b(:) = 0._dp
1067 IF (use_virial) THEN
1068 my_virial_a = 0.0_dp
1069 my_virial_b = 0.0_dp
1070 END IF
1071 END IF
1072 m1 = maxval(npgfa(1:nseta))
1073 ALLOCATE (map_it2(m1, m1))
1074 DO iset = 1, nseta
1075 sgfa = first_sgfa(1, iset)
1076 ncoa = npgfa(iset)*ncoset(la_max(iset))
1077 DO jset = 1, nseta
1078 sgfb = first_sgfa(1, jset)
1079 ncob = npgfa(jset)*ncoset(la_max(jset))
1080 !
1081 map_it2 = .false.
1082 DO ipgf = 1, npgfa(iset)
1083 DO jpgf = 1, npgfa(jset)
1084 zetp = zeta(ipgf, iset) + zeta(jpgf, jset)
1085 igrid_level = gaussian_gridlevel(gridlevel_info, zetp)
1086 rs_grid => rs_v(igrid_level)
1087 map_it2(ipgf, jpgf) = map_gaussian_here(rs_grid, cell%h_inv, ra, offset, group_size, my_pos)
1088 END DO
1089 END DO
1090 offset = offset + 1
1091 !
1092 IF (any(map_it2(1:npgfa(iset), 1:npgfa(jset)))) THEN
1093 hab(:, :) = 0._dp
1094 IF (calculate_forces) THEN
1095 CALL dgemm("N", "N", ncoa, nsgf_seta(jset), nsgf_seta(iset), &
1096 1.0_dp, sphi_a(1, sgfa), SIZE(sphi_a, 1), &
1097 pblk(sgfa, sgfb), SIZE(pblk, 1), &
1098 0.0_dp, work(1, 1), SIZE(work, 1))
1099 CALL dgemm("N", "T", ncoa, ncob, nsgf_seta(jset), &
1100 1.0_dp, work(1, 1), SIZE(work, 1), &
1101 sphi_a(1, sgfb), SIZE(sphi_a, 1), &
1102 0.0_dp, pab(1, 1), SIZE(pab, 1))
1103 END IF
1104
1105 DO ipgf = 1, npgfa(iset)
1106 na1 = (ipgf - 1)*ncoset(la_max(iset))
1107 na2 = ipgf*ncoset(la_max(iset))
1108 DO jpgf = 1, npgfa(jset)
1109 nb1 = (jpgf - 1)*ncoset(la_max(jset))
1110 nb2 = jpgf*ncoset(la_max(jset))
1111 zetp = zeta(ipgf, iset) + zeta(jpgf, jset)
1112 igrid_level = gaussian_gridlevel(gridlevel_info, zetp)
1113 rs_grid => rs_v(igrid_level)
1114
1115 radius = exp_radius_very_extended(la_min=la_min(iset), la_max=la_max(iset), &
1116 lb_min=la_min(jset), lb_max=la_max(jset), &
1117 ra=ra, rb=ra, rp=ra, &
1118 zetp=zetp, eps=eps_rho_rspace, &
1119 prefactor=1.0_dp, cutoff=1.0_dp)
1120
1121 IF (map_it2(ipgf, jpgf)) THEN
1122 IF (calculate_forces) THEN
1123 CALL integrate_pgf_product( &
1124 la_max(iset), zeta(ipgf, iset), la_min(iset), &
1125 la_max(jset), zeta(jpgf, jset), la_min(jset), &
1126 ra=ra, rab=[0.0_dp, 0.0_dp, 0.0_dp], &
1127 rsgrid=rs_v(igrid_level), &
1128 hab=hab, pab=pab, o1=na1, o2=nb1, &
1129 radius=radius, &
1130 calculate_forces=.true., &
1131 force_a=force_a, force_b=force_b, &
1132 use_virial=use_virial, my_virial_a=my_virial_a, my_virial_b=my_virial_b)
1133 ELSE
1134 CALL integrate_pgf_product( &
1135 la_max(iset), zeta(ipgf, iset), la_min(iset), &
1136 la_max(jset), zeta(jpgf, jset), la_min(jset), &
1137 ra=ra, rab=[0.0_dp, 0.0_dp, 0.0_dp], &
1138 rsgrid=rs_v(igrid_level), &
1139 hab=hab, o1=na1, o2=nb1, &
1140 radius=radius, &
1141 calculate_forces=.false.)
1142 END IF
1143 END IF
1144 END DO
1145 END DO
1146 ! contract hab
1147 CALL dgemm("N", "N", ncoa, nsgf_seta(jset), ncob, &
1148 1.0_dp, hab(1, 1), SIZE(hab, 1), &
1149 sphi_a(1, sgfb), SIZE(sphi_a, 1), &
1150 0.0_dp, work(1, 1), SIZE(work, 1))
1151 CALL dgemm("T", "N", nsgf_seta(iset), nsgf_seta(jset), ncoa, &
1152 1.0_dp, sphi_a(1, sgfa), SIZE(sphi_a, 1), &
1153 work(1, 1), SIZE(work, 1), &
1154 1.0_dp, hmat(sgfa, sgfb), SIZE(hmat, 1))
1155 END IF
1156 END DO
1157 END DO
1158 DEALLOCATE (map_it2)
1159 ! update KS matrix
1160 CALL para_env%sum(hmat)
1161 CALL dbcsr_get_block_p(matrix=ksmat, row=atom_a, col=atom_a, block=h_block, found=found)
1162 IF (found) THEN
1163 h_block(1:nsgfa, 1:nsgfa) = h_block(1:nsgfa, 1:nsgfa) + hmat(1:nsgfa, 1:nsgfa)
1164 END IF
1165 IF (calculate_forces) THEN
1166 force(ikind)%rho_elec(:, iatom) = force(ikind)%rho_elec(:, iatom) + 2.0_dp*force_a(:)
1167 IF (use_virial) THEN
1168 IF (use_virial .AND. calculate_forces) THEN
1169 virial%pv_lrigpw = virial%pv_lrigpw + 2.0_dp*my_virial_a
1170 virial%pv_virial = virial%pv_virial + 2.0_dp*my_virial_a
1171 END IF
1172 END IF
1173 END IF
1174 END DO
1175 DEALLOCATE (hab, work, hmat)
1176 IF (calculate_forces) DEALLOCATE (pab, pblk)
1177 END DO
1178
1179 CALL timestop(handle)
1180
1181 END SUBROUTINE integrate_v_rspace_diagonal
1182
1183! **************************************************************************************************
1184!> \brief computes integrals of product of v_rspace times a basis function (vector function)
1185!> and possible forces
1186!> \param qs_env ...
1187!> \param v_rspace ...
1188!> \param f_coef ...
1189!> \param f_integral ...
1190!> \param calculate_forces ...
1191!> \param basis_type ...
1192!> \author JGH [8.2024]
1193! **************************************************************************************************
1194 SUBROUTINE integrate_function(qs_env, v_rspace, f_coef, f_integral, calculate_forces, basis_type)
1195 TYPE(qs_environment_type), POINTER :: qs_env
1196 TYPE(pw_r3d_rs_type), INTENT(IN) :: v_rspace
1197 REAL(kind=dp), DIMENSION(:), INTENT(IN) :: f_coef
1198 REAL(kind=dp), DIMENSION(:), INTENT(OUT) :: f_integral
1199 LOGICAL, INTENT(IN) :: calculate_forces
1200 CHARACTER(len=*), INTENT(IN) :: basis_type
1201
1202 CHARACTER(len=*), PARAMETER :: routinen = 'integrate_function'
1203
1204 INTEGER :: atom_a, group_size, handle, i, iatom, igrid_level, ikind, ipgf, iset, maxco, &
1205 maxsgf_set, my_pos, na1, natom, ncoa, nkind, nseta, offset, sgfa
1206 INTEGER, ALLOCATABLE, DIMENSION(:) :: atom_of_kind
1207 INTEGER, DIMENSION(:), POINTER :: la_max, la_min, npgfa, nsgfa
1208 INTEGER, DIMENSION(:, :), POINTER :: first_sgfa
1209 LOGICAL :: use_virial
1210 REAL(kind=dp) :: eps_rho_rspace, radius
1211 REAL(kind=dp), DIMENSION(3) :: force_a, force_b, ra
1212 REAL(kind=dp), DIMENSION(3, 3) :: my_virial_a, my_virial_b
1213 REAL(kind=dp), DIMENSION(:, :), POINTER :: hab, pab, sphi_a, work, zeta
1214 TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
1215 TYPE(cell_type), POINTER :: cell
1216 TYPE(dft_control_type), POINTER :: dft_control
1217 TYPE(gridlevel_info_type), POINTER :: gridlevel_info
1218 TYPE(gto_basis_set_type), POINTER :: orb_basis_set
1219 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
1220 TYPE(pw_env_type), POINTER :: pw_env
1221 TYPE(qs_force_type), DIMENSION(:), POINTER :: force
1222 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
1223 TYPE(realspace_grid_type), DIMENSION(:), POINTER :: rs_v
1224 TYPE(realspace_grid_type), POINTER :: rs_grid
1225 TYPE(virial_type), POINTER :: virial
1226
1227 CALL timeset(routinen, handle)
1228
1229 CALL get_qs_env(qs_env=qs_env, pw_env=pw_env)
1230 gridlevel_info => pw_env%gridlevel_info
1231
1232 CALL pw_env_get(pw_env, rs_grids=rs_v)
1233 CALL potential_pw2rs(rs_v, v_rspace, pw_env)
1234
1235 CALL get_qs_env(qs_env=qs_env, &
1236 atomic_kind_set=atomic_kind_set, &
1237 qs_kind_set=qs_kind_set, &
1238 force=force, &
1239 cell=cell, &
1240 dft_control=dft_control, &
1241 nkind=nkind, &
1242 natom=natom, &
1243 particle_set=particle_set, &
1244 pw_env=pw_env, &
1245 virial=virial)
1246 CALL get_atomic_kind_set(atomic_kind_set=atomic_kind_set, atom_of_kind=atom_of_kind)
1247
1248 use_virial = virial%pv_availability .AND. (.NOT. virial%pv_numer)
1249 IF (use_virial) THEN
1250 cpabort("Virial NYA")
1251 END IF
1252
1253 eps_rho_rspace = dft_control%qs_control%eps_rho_rspace
1254
1255 CALL get_qs_kind_set(qs_kind_set, &
1256 maxco=maxco, maxsgf_set=maxsgf_set, basis_type=basis_type)
1257 ALLOCATE (hab(maxco, 1), pab(maxco, 1), work(maxco, 1))
1258
1259 offset = 0
1260 my_pos = v_rspace%pw_grid%para%group%mepos
1261 group_size = v_rspace%pw_grid%para%group%num_pe
1262
1263 DO iatom = 1, natom
1264 ikind = particle_set(iatom)%atomic_kind%kind_number
1265 atom_a = atom_of_kind(iatom)
1266 CALL get_qs_kind(qs_kind_set(ikind), basis_set=orb_basis_set, basis_type=basis_type)
1267 CALL get_gto_basis_set(gto_basis_set=orb_basis_set, &
1268 first_sgf=first_sgfa, &
1269 lmax=la_max, &
1270 lmin=la_min, &
1271 npgf=npgfa, &
1272 nset=nseta, &
1273 nsgf_set=nsgfa, &
1274 sphi=sphi_a, &
1275 zet=zeta)
1276 ra(:) = pbc(particle_set(iatom)%r, cell)
1277
1278 force_a(:) = 0._dp
1279 force_b(:) = 0._dp
1280 my_virial_a(:, :) = 0._dp
1281 my_virial_b(:, :) = 0._dp
1282
1283 DO iset = 1, nseta
1284
1285 ncoa = npgfa(iset)*ncoset(la_max(iset))
1286 sgfa = first_sgfa(1, iset)
1287
1288 hab = 0._dp
1289 pab = 0._dp
1290
1291 DO i = 1, nsgfa(iset)
1292 work(i, 1) = f_coef(offset + i)
1293 END DO
1294
1295 CALL dgemm("N", "N", ncoa, 1, nsgfa(iset), &
1296 1.0_dp, sphi_a(1, sgfa), SIZE(sphi_a, 1), &
1297 work(1, 1), SIZE(work, 1), &
1298 0.0_dp, pab(1, 1), SIZE(pab, 1))
1299
1300 DO ipgf = 1, npgfa(iset)
1301
1302 na1 = (ipgf - 1)*ncoset(la_max(iset))
1303
1304 igrid_level = gaussian_gridlevel(gridlevel_info, zeta(ipgf, iset))
1305 rs_grid => rs_v(igrid_level)
1306
1307 IF (map_gaussian_here(rs_grid, cell%h_inv, ra, offset, group_size, my_pos)) THEN
1308 radius = exp_radius_very_extended(la_min=la_min(iset), la_max=la_max(iset), &
1309 lb_min=0, lb_max=0, ra=ra, rb=ra, rp=ra, &
1310 zetp=zeta(ipgf, iset), eps=eps_rho_rspace, &
1311 prefactor=1.0_dp, cutoff=1.0_dp)
1312
1313 IF (calculate_forces) THEN
1314 CALL integrate_pgf_product(la_max=la_max(iset), &
1315 zeta=zeta(ipgf, iset), la_min=la_min(iset), &
1316 lb_max=0, zetb=0.0_dp, lb_min=0, &
1317 ra=ra, rab=[0.0_dp, 0.0_dp, 0.0_dp], &
1318 rsgrid=rs_grid, &
1319 hab=hab, pab=pab, o1=na1, o2=0, radius=radius, &
1320 calculate_forces=.true., &
1321 force_a=force_a, force_b=force_b, &
1322 use_virial=use_virial, &
1323 my_virial_a=my_virial_a, my_virial_b=my_virial_b)
1324 ELSE
1325 CALL integrate_pgf_product(la_max=la_max(iset), &
1326 zeta=zeta(ipgf, iset), la_min=la_min(iset), &
1327 lb_max=0, zetb=0.0_dp, lb_min=0, &
1328 ra=ra, rab=[0.0_dp, 0.0_dp, 0.0_dp], &
1329 rsgrid=rs_grid, &
1330 hab=hab, o1=na1, o2=0, radius=radius, &
1331 calculate_forces=.false.)
1332 END IF
1333
1334 END IF
1335
1336 END DO
1337 !
1338 CALL dgemm("T", "N", nsgfa(iset), 1, ncoa, 1.0_dp, sphi_a(1, sgfa), &
1339 SIZE(sphi_a, 1), hab(1, 1), SIZE(hab, 1), 0.0_dp, work(1, 1), SIZE(work, 1))
1340 DO i = 1, nsgfa(iset)
1341 f_integral(offset + i) = work(i, 1)
1342 END DO
1343
1344 offset = offset + nsgfa(iset)
1345
1346 END DO
1347
1348 IF (calculate_forces) THEN
1349 force(ikind)%rho_elec(:, atom_a) = force(ikind)%rho_elec(:, atom_a) + force_a(:)
1350 IF (use_virial) THEN
1351 virial%pv_virial = virial%pv_virial + my_virial_a
1352 END IF
1353 END IF
1354
1355 END DO
1356
1357 DEALLOCATE (hab, pab, work)
1358
1359 CALL timestop(handle)
1360
1361 END SUBROUTINE integrate_function
1362
1363END MODULE qs_integrate_potential_single
modulo
static GRID_HOST_DEVICE int modulo(int a, int m)
Equivalent of Fortran's MODULO, which always return a positive number. https://gcc....
Definition grid_common.h:125
dgemm
static void dgemm(const char transa, const char transb, const int m, const int n, const int k, const double alpha, const double *a, const int lda, const double *b, const int ldb, const double beta, double *c, const int ldc)
Convenient wrapper to hide Fortran nature of dgemm_, swapping a and b.
Definition grid_cpu_task_list.c:225
cell_types::pbc
Definition cell_types.F:103
external_potential_types::get_potential
Definition external_potential_types.F:278
memory_utilities::reallocate
Definition memory_utilities.F:27
ao_util
All kind of helpful little routines.
Definition ao_util.F:14
ao_util::exp_radius_very_extended
real(kind=dp) function, public exp_radius_very_extended(la_min, la_max, lb_min, lb_max, pab, o1, o2, ra, rb, rp, zetp, eps, prefactor, cutoff, epsabs)
computes the radius of the Gaussian outside of which it is smaller than eps
Definition ao_util.F:209
atomic_kind_types
Define the atomic kind types and their sub types.
Definition atomic_kind_types.F:23
atomic_kind_types::get_atomic_kind_set
subroutine, public get_atomic_kind_set(atomic_kind_set, atom_of_kind, kind_of, natom_of_kind, maxatom, natom, nshell, fist_potential_present, shell_present, shell_adiabatic, shell_check_distance, damping_present)
Get attributes of an atomic kind set.
Definition atomic_kind_types.F:223
atomic_kind_types::get_atomic_kind
subroutine, public get_atomic_kind(atomic_kind, fist_potential, element_symbol, name, mass, kind_number, natom, atom_list, rcov, rvdw, z, qeff, apol, cpol, mm_radius, shell, shell_active, damping)
Get attributes of an atomic kind.
Definition atomic_kind_types.F:138
atprop_types
Holds information on atomic properties.
Definition atprop_types.F:14
atprop_types::atprop_array_init
subroutine, public atprop_array_init(atarray, natom)
...
Definition atprop_types.F:87
basis_set_types
Definition basis_set_types.F:15
basis_set_types::get_gto_basis_set
subroutine, public get_gto_basis_set(gto_basis_set, name, aliases, norm_type, kind_radius, ncgf, nset, nsgf, cgf_symbol, sgf_symbol, norm_cgf, set_radius, lmax, lmin, lx, ly, lz, m, ncgf_set, npgf, nsgf_set, nshell, cphi, pgf_radius, sphi, scon, zet, first_cgf, first_sgf, l, last_cgf, last_sgf, n, gcc, maxco, maxl, maxpgf, maxsgf_set, maxshell, maxso, nco_sum, npgf_sum, nshell_sum, maxder, short_kind_radius, npgf_seg_sum, ccon)
...
Definition basis_set_types.F:643
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_dbcsr_api::dbcsr_get_block_p
subroutine, public dbcsr_get_block_p(matrix, row, col, block, found, row_size, col_size)
...
Definition cp_dbcsr_api.F:693
external_potential_types
Definition of the atomic potential types.
Definition external_potential_types.F:14
gaussian_gridlevels
Definition gaussian_gridlevels.F:13
gaussian_gridlevels::gaussian_gridlevel
integer function, public gaussian_gridlevel(gridlevel_info, exponent)
...
Definition gaussian_gridlevels.F:148
grid_api
Fortran API for the grid package, which is written in C.
Definition grid_api.F:12
grid_api::integrate_pgf_product
subroutine, public integrate_pgf_product(la_max, zeta, la_min, lb_max, zetb, lb_min, ra, rab, rsgrid, hab, pab, o1, o2, radius, calculate_forces, force_a, force_b, compute_tau, use_virial, my_virial_a, my_virial_b, hdab, hadb, a_hdab, use_subpatch, subpatch_pattern)
low level function to compute matrix elements of primitive gaussian functions
Definition grid_api.F:276
kinds
Defines the basic variable types.
Definition kinds.F:23
kinds::dp
integer, parameter, public dp
Definition kinds.F:34
lri_environment_types
contains the types and subroutines for dealing with the lri_env lri : local resolution of the identit...
Definition lri_environment_types.F:18
memory_utilities
Utility routines for the memory handling.
Definition memory_utilities.F:14
message_passing
Interface to the message passing library MPI.
Definition message_passing.F:23
orbital_pointers
Provides Cartesian and spherical orbital pointers and indices.
Definition orbital_pointers.F:15
orbital_pointers::ncoset
integer, dimension(:), allocatable, public ncoset
Definition orbital_pointers.F:42
orbital_pointers::coset
integer, dimension(:, :, :), allocatable, public coset
Definition orbital_pointers.F:45
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_types
Definition pw_types.F:24
qs_environment_types
Definition qs_environment_types.F:14
qs_environment_types::get_qs_env
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, mimic, 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, xcint_weights, 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.
Definition qs_environment_types.F:533
qs_force_types
Definition qs_force_types.F:13
qs_integrate_potential_single
Build up the plane wave density by collocating the primitive Gaussian functions (pgf).
Definition qs_integrate_potential_single.F:30
qs_integrate_potential_single::integrate_function
subroutine, public integrate_function(qs_env, v_rspace, f_coef, f_integral, calculate_forces, basis_type)
computes integrals of product of v_rspace times a basis function (vector function) and possible force...
Definition qs_integrate_potential_single.F:1195
qs_integrate_potential_single::integrate_v_gaussian_rspace
subroutine, public integrate_v_gaussian_rspace(v_rspace, qs_env, alpha, ccore, atecc)
computes the overlap of a set of Gaussians with a potential on grid
Definition qs_integrate_potential_single.F:653
qs_integrate_potential_single::integrate_v_rspace_diagonal
subroutine, public integrate_v_rspace_diagonal(v_rspace, ksmat, pmat, qs_env, calculate_forces, basis_type)
computes integrals of product of v_rspace times the diagonal block basis functions required for LRIGP...
Definition qs_integrate_potential_single.F:978
qs_integrate_potential_single::integrate_v_core_rspace
subroutine, public integrate_v_core_rspace(v_rspace, qs_env, atecc)
computes the forces/virial due to the ionic cores with a potential on grid
Definition qs_integrate_potential_single.F:490
qs_integrate_potential_single::integrate_v_rspace_one_center
subroutine, public integrate_v_rspace_one_center(v_rspace, qs_env, int_res, calculate_forces, basis_type, atomlist)
computes integrals of product of v_rspace times a one-center function required for LRIGPW
Definition qs_integrate_potential_single.F:771
qs_integrate_potential_single::integrate_ppl_rspace
subroutine, public integrate_ppl_rspace(rho_rspace, qs_env)
computes the forces/virial due to the local pseudopotential
Definition qs_integrate_potential_single.F:101
qs_integrate_potential_single::integrate_rho_nlcc
subroutine, public integrate_rho_nlcc(rho_rspace, qs_env)
computes the forces/virial due to the nlcc pseudopotential
Definition qs_integrate_potential_single.F:288
qs_kind_types
Define the quickstep kind type and their sub types.
Definition qs_kind_types.F:23
qs_kind_types::get_qs_kind
subroutine, public get_qs_kind(qs_kind, basis_set, basis_type, ncgf, nsgf, all_potential, tnadd_potential, gth_potential, sgp_potential, upf_potential, cneo_potential, se_parameter, dftb_parameter, xtb_parameter, dftb3_param, zatom, zeff, elec_conf, mao, lmax_dftb, alpha_core_charge, ccore_charge, core_charge, core_charge_radius, paw_proj_set, paw_atom, hard_radius, hard0_radius, max_rad_local, covalent_radius, vdw_radius, gpw_type_forced, harmonics, max_iso_not0, max_s_harm, grid_atom, ngrid_ang, ngrid_rad, lmax_rho0, dft_plus_u_atom, l_of_dft_plus_u, n_of_dft_plus_u, u_minus_j, u_of_dft_plus_u, j_of_dft_plus_u, alpha_of_dft_plus_u, beta_of_dft_plus_u, j0_of_dft_plus_u, occupation_of_dft_plus_u, dispersion, bs_occupation, magnetization, no_optimize, addel, laddel, naddel, orbitals, max_scf, eps_scf, smear, u_ramping, u_minus_j_target, eps_u_ramping, init_u_ramping_each_scf, reltmat, ghost, monovalent, floating, name, element_symbol, pao_basis_size, pao_model_file, pao_potentials, pao_descriptors, nelec)
Get attributes of an atomic kind.
Definition qs_kind_types.F:468
qs_kind_types::get_qs_kind_set
subroutine, public get_qs_kind_set(qs_kind_set, all_potential_present, tnadd_potential_present, gth_potential_present, sgp_potential_present, paw_atom_present, dft_plus_u_atom_present, maxcgf, maxsgf, maxco, maxco_proj, maxgtops, maxlgto, maxlprj, maxnset, maxsgf_set, ncgf, npgf, nset, nsgf, nshell, maxpol, maxlppl, maxlppnl, maxppnl, nelectron, maxder, max_ngrid_rad, max_sph_harm, maxg_iso_not0, lmax_rho0, basis_rcut, basis_type, total_zeff_corr, npgf_seg, cneo_potential_present, nkind_q, natom_q)
Get attributes of an atomic kind set.
Definition qs_kind_types.F:913
realspace_grid_types
Definition realspace_grid_types.F:18
realspace_grid_types::transfer_pw2rs
subroutine, public transfer_pw2rs(rs, pw)
...
Definition realspace_grid_types.F:695
realspace_grid_types::map_gaussian_here
pure logical function, public map_gaussian_here(rs_grid, h_inv, ra, offset, group_size, my_pos)
...
Definition realspace_grid_types.F:2146
realspace_grid_types::rs_grid_zero
subroutine, public rs_grid_zero(rs)
Initialize grid to zero.
Definition realspace_grid_types.F:1907
rs_pw_interface
Transfers densities from PW to RS grids and potentials from PW to RS.
Definition rs_pw_interface.F:14
rs_pw_interface::potential_pw2rs
subroutine, public potential_pw2rs(rs_v, v_rspace, pw_env)
transfers a potential from a pw_grid to a vector of realspace multigrids
Definition rs_pw_interface.F:152
virial_types
Definition virial_types.F:13
atomic_kind_types::atomic_kind_type
Provides all information about an atomic kind.
Definition atomic_kind_types.F:45
atprop_types::atprop_type
type for the atomic properties
Definition atprop_types.F:31
basis_set_types::gto_basis_set_type
Definition basis_set_types.F:72
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:757
cp_dbcsr_api::dbcsr_type
Definition cp_dbcsr_api.F:176
external_potential_types::gth_potential_type
Definition external_potential_types.F:116
gaussian_gridlevels::gridlevel_info_type
Definition gaussian_gridlevels.F:34
lri_environment_types::lri_kind_type
Definition lri_environment_types.F:332
message_passing::mp_para_env_type
stores all the informations relevant to an mpi environment
Definition message_passing.F:746
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_types::pw_r3d_rs_type
Definition pw_types.F:62
qs_environment_types::qs_environment_type
Definition qs_environment_types.F:222
qs_force_types::qs_force_type
Definition qs_force_types.F:28
qs_kind_types::qs_kind_type
Provides all information about a quickstep kind.
Definition qs_kind_types.F:177
realspace_grid_types::realspace_grid_type
Definition realspace_grid_types.F:136
virial_types::virial_type
Definition virial_types.F:26