(git:2d33723)
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 cpabort("")
214 END SELECT
215 END DO
216
217 DO iatom = 1, natom_of_kind
218 atom_a = atom_list(iatom)
219 ra(:) = pbc(particle_set(atom_a)%r, cell)
220 IF (rs_v%desc%parallel .AND. .NOT. rs_v%desc%distributed) THEN
221 ! replicated realspace grid, split the atoms up between procs
222 IF (modulo(iatom, rs_v%desc%group_size) == rs_v%desc%my_pos) THEN
223 npme = npme + 1
224 cores(npme) = iatom
225 END IF
226 ELSE
227 npme = npme + 1
228 cores(npme) = iatom
229 END IF
230 END DO
231
232 DO j = 1, npme
233
234 iatom = cores(j)
235 atom_a = atom_list(iatom)
236 ra(:) = pbc(particle_set(atom_a)%r, cell)
237 hab(:, 1) = 0.0_dp
238 force_a(:) = 0.0_dp
239 force_b(:) = 0.0_dp
240 IF (use_virial) THEN
241 my_virial_a = 0.0_dp
242 my_virial_b = 0.0_dp
243 END IF
244 ni = 2*lppl - 2
245
246 radius = exp_radius_very_extended(la_min=0, la_max=ni, lb_min=0, lb_max=0, &
247 ra=ra, rb=ra, rp=ra, &
248 zetp=alpha, eps=eps_rho_rspace, &
249 pab=pab, o1=0, o2=0, & ! without map_consistent
250 prefactor=1.0_dp, cutoff=1.0_dp)
251
252 CALL integrate_pgf_product(ni, alpha, 0, &
253 0, 0.0_dp, 0, ra, [0.0_dp, 0.0_dp, 0.0_dp], &
254 rs_v, hab, pab=pab, o1=0, o2=0, &
255 radius=radius, &
256 calculate_forces=.true., force_a=force_a, &
257 force_b=force_b, use_virial=use_virial, my_virial_a=my_virial_a, &
258 my_virial_b=my_virial_b, use_subpatch=.true., subpatch_pattern=0)
259
260 force(ikind)%gth_ppl(:, iatom) = &
261 force(ikind)%gth_ppl(:, iatom) + force_a(:)*rho_rspace%pw_grid%dvol
262
263 IF (use_virial) THEN
264 virial%pv_ppl = virial%pv_ppl + my_virial_a*rho_rspace%pw_grid%dvol
265 virial%pv_virial = virial%pv_virial + my_virial_a*rho_rspace%pw_grid%dvol
266 cpabort("Virial not debuged for CORE_PPL")
267 END IF
268 END DO
269
270 DEALLOCATE (hab, pab)
271
272 END DO
273
274 DEALLOCATE (cores)
275
276 CALL timestop(handle)
277
278 END SUBROUTINE integrate_ppl_rspace
279
280! **************************************************************************************************
281!> \brief computes the forces/virial due to the nlcc pseudopotential
282!> \param rho_rspace ...
283!> \param qs_env ...
284! **************************************************************************************************
285 SUBROUTINE integrate_rho_nlcc(rho_rspace, qs_env)
286 TYPE(pw_r3d_rs_type), INTENT(IN) :: rho_rspace
287 TYPE(qs_environment_type), POINTER :: qs_env
288
289 CHARACTER(len=*), PARAMETER :: routinen = 'integrate_rho_nlcc'
290
291 INTEGER :: atom_a, handle, iatom, iexp_nlcc, ikind, &
292 ithread, j, n, natom, nc, nexp_nlcc, &
293 ni, npme, nthread
294 INTEGER, DIMENSION(:), POINTER :: atom_list, cores, nct_nlcc
295 LOGICAL :: nlcc, use_virial
296 REAL(kind=dp) :: alpha, eps_rho_rspace, radius
297 REAL(kind=dp), DIMENSION(3) :: force_a, force_b, ra
298 REAL(kind=dp), DIMENSION(3, 3) :: my_virial_a, my_virial_b
299 REAL(kind=dp), DIMENSION(:), POINTER :: alpha_nlcc
300 REAL(kind=dp), DIMENSION(:, :), POINTER :: cval_nlcc, hab, pab
301 TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
302 TYPE(cell_type), POINTER :: cell
303 TYPE(dft_control_type), POINTER :: dft_control
304 TYPE(gth_potential_type), POINTER :: gth_potential
305 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
306 TYPE(pw_env_type), POINTER :: pw_env
307 TYPE(qs_force_type), DIMENSION(:), POINTER :: force
308 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
309 TYPE(realspace_grid_type), POINTER :: rs_v
310 TYPE(virial_type), POINTER :: virial
311
312 CALL timeset(routinen, handle)
313
314 NULLIFY (pw_env, cores)
315
316 CALL get_qs_env(qs_env=qs_env, pw_env=pw_env)
317 CALL pw_env_get(pw_env=pw_env, auxbas_rs_grid=rs_v)
318
319 CALL transfer_pw2rs(rs_v, rho_rspace)
320
321 CALL get_qs_env(qs_env=qs_env, &
322 atomic_kind_set=atomic_kind_set, &
323 qs_kind_set=qs_kind_set, &
324 cell=cell, &
325 dft_control=dft_control, &
326 particle_set=particle_set, &
327 pw_env=pw_env, &
328 force=force, virial=virial)
329
330 use_virial = virial%pv_availability .AND. (.NOT. virial%pv_numer)
331
332 eps_rho_rspace = dft_control%qs_control%eps_rho_rspace
333
334 DO ikind = 1, SIZE(atomic_kind_set)
335
336 CALL get_atomic_kind(atomic_kind_set(ikind), natom=natom, atom_list=atom_list)
337 CALL get_qs_kind(qs_kind_set(ikind), gth_potential=gth_potential)
338
339 IF (.NOT. ASSOCIATED(gth_potential)) cycle
340 CALL get_potential(potential=gth_potential, nlcc_present=nlcc, nexp_nlcc=nexp_nlcc, &
341 alpha_nlcc=alpha_nlcc, nct_nlcc=nct_nlcc, cval_nlcc=cval_nlcc)
342
343 IF (.NOT. nlcc) cycle
344
345 DO iexp_nlcc = 1, nexp_nlcc
346
347 alpha = alpha_nlcc(iexp_nlcc)
348 nc = nct_nlcc(iexp_nlcc)
349
350 ni = ncoset(2*nc - 2)
351
352 nthread = 1
353 ithread = 0
354
355 ALLOCATE (hab(ni, 1), pab(ni, 1))
356 pab = 0._dp
357
358 CALL reallocate(cores, 1, natom)
359 npme = 0
360 cores = 0
361
362 ! prepare core function
363 DO j = 1, nc
364 SELECT CASE (j)
365 CASE (1)
366 pab(1, 1) = cval_nlcc(1, iexp_nlcc)
367 CASE (2)
368 n = coset(2, 0, 0)
369 pab(n, 1) = cval_nlcc(2, iexp_nlcc)/alpha**2
370 n = coset(0, 2, 0)
371 pab(n, 1) = cval_nlcc(2, iexp_nlcc)/alpha**2
372 n = coset(0, 0, 2)
373 pab(n, 1) = cval_nlcc(2, iexp_nlcc)/alpha**2
374 CASE (3)
375 n = coset(4, 0, 0)
376 pab(n, 1) = cval_nlcc(3, iexp_nlcc)/alpha**4
377 n = coset(0, 4, 0)
378 pab(n, 1) = cval_nlcc(3, iexp_nlcc)/alpha**4
379 n = coset(0, 0, 4)
380 pab(n, 1) = cval_nlcc(3, iexp_nlcc)/alpha**4
381 n = coset(2, 2, 0)
382 pab(n, 1) = 2._dp*cval_nlcc(3, iexp_nlcc)/alpha**4
383 n = coset(2, 0, 2)
384 pab(n, 1) = 2._dp*cval_nlcc(3, iexp_nlcc)/alpha**4
385 n = coset(0, 2, 2)
386 pab(n, 1) = 2._dp*cval_nlcc(3, iexp_nlcc)/alpha**4
387 CASE (4)
388 n = coset(6, 0, 0)
389 pab(n, 1) = cval_nlcc(4, iexp_nlcc)/alpha**6
390 n = coset(0, 6, 0)
391 pab(n, 1) = cval_nlcc(4, iexp_nlcc)/alpha**6
392 n = coset(0, 0, 6)
393 pab(n, 1) = cval_nlcc(4, iexp_nlcc)/alpha**6
394 n = coset(4, 2, 0)
395 pab(n, 1) = 3._dp*cval_nlcc(4, iexp_nlcc)/alpha**6
396 n = coset(4, 0, 2)
397 pab(n, 1) = 3._dp*cval_nlcc(4, iexp_nlcc)/alpha**6
398 n = coset(2, 4, 0)
399 pab(n, 1) = 3._dp*cval_nlcc(4, iexp_nlcc)/alpha**6
400 n = coset(2, 0, 4)
401 pab(n, 1) = 3._dp*cval_nlcc(4, iexp_nlcc)/alpha**6
402 n = coset(0, 4, 2)
403 pab(n, 1) = 3._dp*cval_nlcc(4, iexp_nlcc)/alpha**6
404 n = coset(0, 2, 4)
405 pab(n, 1) = 3._dp*cval_nlcc(4, iexp_nlcc)/alpha**6
406 n = coset(2, 2, 2)
407 pab(n, 1) = 6._dp*cval_nlcc(4, iexp_nlcc)/alpha**6
408 CASE DEFAULT
409 cpabort("")
410 END SELECT
411 END DO
412 IF (dft_control%nspins == 2) pab = pab*0.5_dp
413
414 DO iatom = 1, natom
415 atom_a = atom_list(iatom)
416 ra(:) = pbc(particle_set(atom_a)%r, cell)
417 IF (rs_v%desc%parallel .AND. .NOT. rs_v%desc%distributed) THEN
418 ! replicated realspace grid, split the atoms up between procs
419 IF (modulo(iatom, rs_v%desc%group_size) == rs_v%desc%my_pos) THEN
420 npme = npme + 1
421 cores(npme) = iatom
422 END IF
423 ELSE
424 npme = npme + 1
425 cores(npme) = iatom
426 END IF
427 END DO
428
429 DO j = 1, npme
430
431 iatom = cores(j)
432 atom_a = atom_list(iatom)
433 ra(:) = pbc(particle_set(atom_a)%r, cell)
434 hab(:, 1) = 0.0_dp
435 force_a(:) = 0.0_dp
436 force_b(:) = 0.0_dp
437 IF (use_virial) THEN
438 my_virial_a = 0.0_dp
439 my_virial_b = 0.0_dp
440 END IF
441 ni = 2*nc - 2
442
443 radius = exp_radius_very_extended(la_min=0, la_max=ni, lb_min=0, lb_max=0, &
444 ra=ra, rb=ra, rp=ra, &
445 zetp=1/(2*alpha**2), eps=eps_rho_rspace, &
446 pab=pab, o1=0, o2=0, & ! without map_consistent
447 prefactor=1.0_dp, cutoff=1.0_dp)
448
449 CALL integrate_pgf_product(ni, 1/(2*alpha**2), 0, &
450 0, 0.0_dp, 0, ra, [0.0_dp, 0.0_dp, 0.0_dp], &
451 rs_v, hab, pab=pab, o1=0, o2=0, &
452 radius=radius, &
453 calculate_forces=.true., force_a=force_a, &
454 force_b=force_b, use_virial=use_virial, my_virial_a=my_virial_a, &
455 my_virial_b=my_virial_b, use_subpatch=.true., subpatch_pattern=0)
456
457 force(ikind)%gth_nlcc(:, iatom) = &
458 force(ikind)%gth_nlcc(:, iatom) + force_a(:)*rho_rspace%pw_grid%dvol
459
460 IF (use_virial) THEN
461 virial%pv_nlcc = virial%pv_nlcc + my_virial_a*rho_rspace%pw_grid%dvol
462 virial%pv_virial = virial%pv_virial + my_virial_a*rho_rspace%pw_grid%dvol
463 END IF
464 END DO
465
466 DEALLOCATE (hab, pab)
467
468 END DO
469
470 END DO
471
472 DEALLOCATE (cores)
473
474 CALL timestop(handle)
475
476 END SUBROUTINE integrate_rho_nlcc
477
478! **************************************************************************************************
479!> \brief computes the forces/virial due to the ionic cores with a potential on
480!> grid
481!> \param v_rspace ...
482!> \param qs_env ...
483!> \param atecc ...
484! **************************************************************************************************
485 SUBROUTINE integrate_v_core_rspace(v_rspace, qs_env, atecc)
486 TYPE(pw_r3d_rs_type), INTENT(IN) :: v_rspace
487 TYPE(qs_environment_type), POINTER :: qs_env
488 REAL(kind=dp), DIMENSION(:), OPTIONAL :: atecc
489
490 CHARACTER(len=*), PARAMETER :: routinen = 'integrate_v_core_rspace'
491
492 INTEGER :: atom_a, handle, iatom, ikind, j, natom, &
493 natom_of_kind, npme
494 INTEGER, DIMENSION(:), POINTER :: atom_list, cores
495 LOGICAL :: paw_atom, skip_fcore, use_virial
496 REAL(kind=dp) :: alpha_core_charge, ccore_charge, &
497 eps_rho_rspace, radius
498 REAL(kind=dp), DIMENSION(3) :: force_a, force_b, ra
499 REAL(kind=dp), DIMENSION(3, 3) :: my_virial_a, my_virial_b
500 REAL(kind=dp), DIMENSION(:, :), POINTER :: hab, pab
501 TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
502 TYPE(atprop_type), POINTER :: atprop
503 TYPE(cell_type), POINTER :: cell
504 TYPE(dft_control_type), POINTER :: dft_control
505 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
506 TYPE(pw_env_type), POINTER :: pw_env
507 TYPE(qs_force_type), DIMENSION(:), POINTER :: force
508 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
509 TYPE(realspace_grid_type), POINTER :: rs_v
510 TYPE(virial_type), POINTER :: virial
511
512 CALL timeset(routinen, handle)
513 NULLIFY (virial, force, atprop, dft_control)
514
515 CALL get_qs_env(qs_env=qs_env, dft_control=dft_control)
516
517 !If gapw, check for gpw kinds
518 skip_fcore = .false.
519 IF (dft_control%qs_control%gapw) THEN
520 IF (.NOT. dft_control%qs_control%gapw_control%nopaw_as_gpw) skip_fcore = .true.
521 END IF
522
523 IF (.NOT. skip_fcore) THEN
524 NULLIFY (pw_env)
525 ALLOCATE (cores(1))
526 ALLOCATE (hab(1, 1))
527 ALLOCATE (pab(1, 1))
528
529 CALL get_qs_env(qs_env=qs_env, pw_env=pw_env)
530 CALL pw_env_get(pw_env=pw_env, auxbas_rs_grid=rs_v)
531
532 CALL transfer_pw2rs(rs_v, v_rspace)
533
534 CALL get_qs_env(qs_env=qs_env, &
535 atomic_kind_set=atomic_kind_set, &
536 qs_kind_set=qs_kind_set, &
537 cell=cell, &
538 dft_control=dft_control, &
539 particle_set=particle_set, &
540 pw_env=pw_env, &
541 force=force, &
542 virial=virial, &
543 atprop=atprop)
544
545 ! atomic energy contributions
546 natom = SIZE(particle_set)
547 IF (ASSOCIATED(atprop)) THEN
548 CALL atprop_array_init(atprop%ateb, natom)
549 END IF
550
551 use_virial = virial%pv_availability .AND. (.NOT. virial%pv_numer)
552
553 eps_rho_rspace = dft_control%qs_control%eps_rho_rspace
554
555 DO ikind = 1, SIZE(atomic_kind_set)
556
557 CALL get_atomic_kind(atomic_kind_set(ikind), natom=natom_of_kind, atom_list=atom_list)
558 CALL get_qs_kind(qs_kind_set(ikind), paw_atom=paw_atom)
559
560 IF (dft_control%qs_control%gapw .AND. paw_atom) cycle
561
562 CALL get_qs_kind(qs_kind_set(ikind), alpha_core_charge=alpha_core_charge, &
563 ccore_charge=ccore_charge)
564
565 pab(1, 1) = -ccore_charge
566 IF (alpha_core_charge == 0.0_dp .OR. pab(1, 1) == 0.0_dp) cycle
567
568 CALL reallocate(cores, 1, natom_of_kind)
569 npme = 0
570 cores = 0
571
572 DO iatom = 1, natom_of_kind
573 atom_a = atom_list(iatom)
574 ra(:) = pbc(particle_set(atom_a)%r, cell)
575 IF (rs_v%desc%parallel .AND. .NOT. rs_v%desc%distributed) THEN
576 ! replicated realspace grid, split the atoms up between procs
577 IF (modulo(iatom, rs_v%desc%group_size) == rs_v%desc%my_pos) THEN
578 npme = npme + 1
579 cores(npme) = iatom
580 END IF
581 ELSE
582 npme = npme + 1
583 cores(npme) = iatom
584 END IF
585 END DO
586
587 DO j = 1, npme
588
589 iatom = cores(j)
590 atom_a = atom_list(iatom)
591 ra(:) = pbc(particle_set(atom_a)%r, cell)
592 hab(1, 1) = 0.0_dp
593 force_a(:) = 0.0_dp
594 force_b(:) = 0.0_dp
595 IF (use_virial) THEN
596 my_virial_a = 0.0_dp
597 my_virial_b = 0.0_dp
598 END IF
599
600 radius = exp_radius_very_extended(la_min=0, la_max=0, lb_min=0, lb_max=0, &
601 ra=ra, rb=ra, rp=ra, &
602 zetp=alpha_core_charge, eps=eps_rho_rspace, &
603 pab=pab, o1=0, o2=0, & ! without map_consistent
604 prefactor=1.0_dp, cutoff=1.0_dp)
605
606 CALL integrate_pgf_product(0, alpha_core_charge, 0, &
607 0, 0.0_dp, 0, ra, [0.0_dp, 0.0_dp, 0.0_dp], &
608 rs_v, hab, pab=pab, o1=0, o2=0, &
609 radius=radius, &
610 calculate_forces=.true., force_a=force_a, &
611 force_b=force_b, use_virial=use_virial, my_virial_a=my_virial_a, &
612 my_virial_b=my_virial_b, use_subpatch=.true., subpatch_pattern=0)
613
614 IF (ASSOCIATED(force)) THEN
615 force(ikind)%rho_core(:, iatom) = force(ikind)%rho_core(:, iatom) + force_a(:)
616 END IF
617 IF (use_virial) THEN
618 virial%pv_ehartree = virial%pv_ehartree + my_virial_a
619 virial%pv_virial = virial%pv_virial + my_virial_a
620 END IF
621 IF (ASSOCIATED(atprop)) THEN
622 atprop%ateb(atom_a) = atprop%ateb(atom_a) + 0.5_dp*hab(1, 1)*pab(1, 1)
623 END IF
624 IF (PRESENT(atecc)) THEN
625 atecc(atom_a) = atecc(atom_a) + 0.5_dp*hab(1, 1)*pab(1, 1)
626 END IF
627
628 END DO
629
630 END DO
631
632 DEALLOCATE (hab, pab, cores)
633
634 END IF
635
636 CALL timestop(handle)
637
638 END SUBROUTINE integrate_v_core_rspace
639
640! **************************************************************************************************
641!> \brief computes the overlap of a set of Gaussians with a potential on grid
642!> \param v_rspace ...
643!> \param qs_env ...
644!> \param alpha ...
645!> \param ccore ...
646!> \param atecc ...
647! **************************************************************************************************
648 SUBROUTINE integrate_v_gaussian_rspace(v_rspace, qs_env, alpha, ccore, atecc)
649 TYPE(pw_r3d_rs_type), INTENT(IN) :: v_rspace
650 TYPE(qs_environment_type), POINTER :: qs_env
651 REAL(kind=dp), DIMENSION(:), INTENT(IN) :: alpha, ccore
652 REAL(kind=dp), DIMENSION(:) :: atecc
653
654 CHARACTER(len=*), PARAMETER :: routinen = 'integrate_v_gaussian_rspace'
655
656 INTEGER :: atom_a, handle, iatom, ikind, j, natom, &
657 natom_of_kind, npme
658 INTEGER, DIMENSION(:), POINTER :: atom_list, cores
659 REAL(kind=dp) :: alpha_core_charge, eps_rho_rspace, radius
660 REAL(kind=dp), DIMENSION(3) :: ra
661 REAL(kind=dp), DIMENSION(:, :), POINTER :: hab, pab
662 TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
663 TYPE(cell_type), POINTER :: cell
664 TYPE(dft_control_type), POINTER :: dft_control
665 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
666 TYPE(pw_env_type), POINTER :: pw_env
667 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
668 TYPE(realspace_grid_type), POINTER :: rs_v
669
670 CALL timeset(routinen, handle)
671
672 CALL get_qs_env(qs_env=qs_env, dft_control=dft_control)
673
674 !If gapw, check for gpw kinds
675 cpassert(.NOT. dft_control%qs_control%gapw)
676
677 NULLIFY (pw_env)
678 ALLOCATE (cores(1))
679 ALLOCATE (hab(1, 1))
680 ALLOCATE (pab(1, 1))
681
682 CALL get_qs_env(qs_env=qs_env, pw_env=pw_env)
683 CALL pw_env_get(pw_env=pw_env, auxbas_rs_grid=rs_v)
684
685 CALL transfer_pw2rs(rs_v, v_rspace)
686
687 CALL get_qs_env(qs_env=qs_env, &
688 atomic_kind_set=atomic_kind_set, &
689 qs_kind_set=qs_kind_set, &
690 cell=cell, &
691 dft_control=dft_control, &
692 particle_set=particle_set, &
693 pw_env=pw_env)
694
695 ! atomic energy contributions
696 natom = SIZE(particle_set)
697 eps_rho_rspace = dft_control%qs_control%eps_rho_rspace
698
699 DO ikind = 1, SIZE(atomic_kind_set)
700
701 CALL get_atomic_kind(atomic_kind_set(ikind), natom=natom_of_kind, atom_list=atom_list)
702 pab(1, 1) = -ccore(ikind)
703 alpha_core_charge = alpha(ikind)
704 IF (alpha_core_charge == 0.0_dp .OR. pab(1, 1) == 0.0_dp) cycle
705
706 CALL reallocate(cores, 1, natom_of_kind)
707 npme = 0
708 cores = 0
709
710 DO iatom = 1, natom_of_kind
711 atom_a = atom_list(iatom)
712 ra(:) = pbc(particle_set(atom_a)%r, cell)
713 IF (rs_v%desc%parallel .AND. .NOT. rs_v%desc%distributed) THEN
714 ! replicated realspace grid, split the atoms up between procs
715 IF (modulo(iatom, rs_v%desc%group_size) == rs_v%desc%my_pos) THEN
716 npme = npme + 1
717 cores(npme) = iatom
718 END IF
719 ELSE
720 npme = npme + 1
721 cores(npme) = iatom
722 END IF
723 END DO
724
725 DO j = 1, npme
726
727 iatom = cores(j)
728 atom_a = atom_list(iatom)
729 ra(:) = pbc(particle_set(atom_a)%r, cell)
730 hab(1, 1) = 0.0_dp
731
732 radius = exp_radius_very_extended(la_min=0, la_max=0, lb_min=0, lb_max=0, &
733 ra=ra, rb=ra, rp=ra, &
734 zetp=alpha_core_charge, eps=eps_rho_rspace, &
735 pab=pab, o1=0, o2=0, & ! without map_consistent
736 prefactor=1.0_dp, cutoff=1.0_dp)
737
738 CALL integrate_pgf_product(0, alpha_core_charge, 0, &
739 0, 0.0_dp, 0, ra, [0.0_dp, 0.0_dp, 0.0_dp], &
740 rs_v, hab, pab=pab, o1=0, o2=0, &
741 radius=radius, calculate_forces=.false., &
742 use_subpatch=.true., subpatch_pattern=0)
743 atecc(atom_a) = atecc(atom_a) + 0.5_dp*hab(1, 1)*pab(1, 1)
744
745 END DO
746
747 END DO
748
749 DEALLOCATE (hab, pab, cores)
750
751 CALL timestop(handle)
752
753 END SUBROUTINE integrate_v_gaussian_rspace
754! **************************************************************************************************
755!> \brief computes integrals of product of v_rspace times a one-center function
756!> required for LRIGPW
757!> \param v_rspace ...
758!> \param qs_env ...
759!> \param int_res ...
760!> \param calculate_forces ...
761!> \param basis_type ...
762!> \param atomlist ...
763!> \author Dorothea Golze
764! **************************************************************************************************
765 SUBROUTINE integrate_v_rspace_one_center(v_rspace, qs_env, int_res, &
766 calculate_forces, basis_type, atomlist)
767 TYPE(pw_r3d_rs_type), INTENT(IN) :: v_rspace
768 TYPE(qs_environment_type), POINTER :: qs_env
769 TYPE(lri_kind_type), DIMENSION(:), POINTER :: int_res
770 LOGICAL, INTENT(IN) :: calculate_forces
771 CHARACTER(len=*), INTENT(IN) :: basis_type
772 INTEGER, DIMENSION(:), OPTIONAL :: atomlist
773
774 CHARACTER(len=*), PARAMETER :: routinen = 'integrate_v_rspace_one_center'
775
776 INTEGER :: atom_a, group_size, handle, i, iatom, igrid_level, ikind, ipgf, iset, m1, &
777 max_npgf, maxco, maxsgf_set, my_pos, na1, natom_of_kind, ncoa, nkind, nseta, offset, sgfa
778 INTEGER, DIMENSION(:), POINTER :: atom_list, la_max, la_min, npgfa, &
779 nsgf_seta
780 INTEGER, DIMENSION(:, :), POINTER :: first_sgfa
781 LOGICAL :: use_virial
782 LOGICAL, ALLOCATABLE, DIMENSION(:) :: map_it
783 REAL(kind=dp) :: eps_rho_rspace, radius
784 REAL(kind=dp), DIMENSION(3) :: force_a, force_b, ra
785 REAL(kind=dp), DIMENSION(3, 3) :: my_virial_a, my_virial_b
786 REAL(kind=dp), DIMENSION(:), POINTER :: set_radius_a
787 REAL(kind=dp), DIMENSION(:, :), POINTER :: hab, pab, rpgfa, sphi_a, work_f, work_i, &
788 zeta
789 TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
790 TYPE(cell_type), POINTER :: cell
791 TYPE(dft_control_type), POINTER :: dft_control
792 TYPE(gridlevel_info_type), POINTER :: gridlevel_info
793 TYPE(gto_basis_set_type), POINTER :: lri_basis_set
794 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
795 TYPE(pw_env_type), POINTER :: pw_env
796 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
797 TYPE(realspace_grid_type), DIMENSION(:), POINTER :: rs_v
798 TYPE(realspace_grid_type), POINTER :: rs_grid
799 TYPE(virial_type), POINTER :: virial
800
801 CALL timeset(routinen, handle)
802
803 NULLIFY (atomic_kind_set, qs_kind_set, atom_list, cell, dft_control, &
804 first_sgfa, gridlevel_info, hab, la_max, la_min, lri_basis_set, &
805 npgfa, nsgf_seta, pab, particle_set, pw_env, rpgfa, &
806 rs_grid, rs_v, virial, set_radius_a, sphi_a, work_f, &
807 work_i, zeta)
808
809 CALL get_qs_env(qs_env=qs_env, pw_env=pw_env)
810
811 CALL pw_env_get(pw_env, rs_grids=rs_v)
812 DO i = 1, SIZE(rs_v)
813 CALL rs_grid_zero(rs_v(i))
814 END DO
815
816 gridlevel_info => pw_env%gridlevel_info
817
818 CALL potential_pw2rs(rs_v, v_rspace, pw_env)
819
820 CALL get_qs_env(qs_env=qs_env, &
821 atomic_kind_set=atomic_kind_set, &
822 qs_kind_set=qs_kind_set, &
823 cell=cell, &
824 dft_control=dft_control, &
825 nkind=nkind, &
826 particle_set=particle_set, &
827 pw_env=pw_env, &
828 virial=virial)
829
830 use_virial = virial%pv_availability .AND. (.NOT. virial%pv_numer)
831
832 eps_rho_rspace = dft_control%qs_control%eps_rho_rspace
833
834 offset = 0
835 my_pos = v_rspace%pw_grid%para%group%mepos
836 group_size = v_rspace%pw_grid%para%group%num_pe
837
838 DO ikind = 1, nkind
839
840 CALL get_atomic_kind(atomic_kind_set(ikind), natom=natom_of_kind, atom_list=atom_list)
841 CALL get_qs_kind(qs_kind_set(ikind), basis_set=lri_basis_set, basis_type=basis_type)
842 CALL get_gto_basis_set(gto_basis_set=lri_basis_set, &
843 first_sgf=first_sgfa, &
844 lmax=la_max, &
845 lmin=la_min, &
846 maxco=maxco, &
847 maxsgf_set=maxsgf_set, &
848 npgf=npgfa, &
849 nset=nseta, &
850 nsgf_set=nsgf_seta, &
851 pgf_radius=rpgfa, &
852 set_radius=set_radius_a, &
853 sphi=sphi_a, &
854 zet=zeta)
855
856 ALLOCATE (hab(maxco, 1), pab(maxco, 1))
857 hab = 0._dp
858 pab(:, 1) = 0._dp
859 max_npgf = maxval(npgfa(1:nseta))
860 ALLOCATE (map_it(max_npgf))
861 ALLOCATE (work_i(maxsgf_set, 1))
862 IF (calculate_forces) ALLOCATE (work_f(maxsgf_set, 1))
863
864 DO iatom = 1, natom_of_kind
865
866 atom_a = atom_list(iatom)
867 IF (PRESENT(atomlist)) THEN
868 IF (atomlist(atom_a) == 0) cycle
869 END IF
870 ra(:) = pbc(particle_set(atom_a)%r, cell)
871 force_a(:) = 0._dp
872 force_b(:) = 0._dp
873 my_virial_a(:, :) = 0._dp
874 my_virial_b(:, :) = 0._dp
875
876 DO iset = 1, nseta
877 !
878 map_it = .false.
879 DO ipgf = 1, npgfa(iset)
880 igrid_level = gaussian_gridlevel(gridlevel_info, zeta(ipgf, iset))
881 rs_grid => rs_v(igrid_level)
882 map_it(ipgf) = map_gaussian_here(rs_grid, cell%h_inv, ra, offset, group_size, my_pos)
883 END DO
884 offset = offset + 1
885 !
886 IF (any(map_it(1:npgfa(iset)))) THEN
887 sgfa = first_sgfa(1, iset)
888 ncoa = npgfa(iset)*ncoset(la_max(iset))
889 hab(:, 1) = 0._dp
890 work_i(1:nsgf_seta(iset), 1) = 0.0_dp
891
892 ! get fit coefficients for forces
893 IF (calculate_forces) THEN
894 m1 = sgfa + nsgf_seta(iset) - 1
895 work_f(1:nsgf_seta(iset), 1) = int_res(ikind)%acoef(iatom, sgfa:m1)
896 CALL dgemm("N", "N", ncoa, 1, nsgf_seta(iset), 1.0_dp, sphi_a(1, sgfa), &
897 SIZE(sphi_a, 1), work_f(1, 1), SIZE(work_f, 1), 0.0_dp, pab(1, 1), &
898 SIZE(pab, 1))
899 END IF
900
901 DO ipgf = 1, npgfa(iset)
902 na1 = (ipgf - 1)*ncoset(la_max(iset))
903 igrid_level = gaussian_gridlevel(gridlevel_info, zeta(ipgf, iset))
904 rs_grid => rs_v(igrid_level)
905
906 radius = exp_radius_very_extended(la_min=la_min(iset), la_max=la_max(iset), &
907 lb_min=0, lb_max=0, ra=ra, rb=ra, rp=ra, &
908 zetp=zeta(ipgf, iset), eps=eps_rho_rspace, &
909 prefactor=1.0_dp, cutoff=1.0_dp)
910
911 IF (map_it(ipgf)) THEN
912 IF (.NOT. calculate_forces) THEN
913 CALL integrate_pgf_product(la_max=la_max(iset), &
914 zeta=zeta(ipgf, iset), la_min=la_min(iset), &
915 lb_max=0, zetb=0.0_dp, lb_min=0, &
916 ra=ra, rab=[0.0_dp, 0.0_dp, 0.0_dp], &
917 rsgrid=rs_grid, &
918 hab=hab, o1=na1, o2=0, radius=radius, &
919 calculate_forces=calculate_forces)
920 ELSE
921 CALL integrate_pgf_product(la_max=la_max(iset), &
922 zeta=zeta(ipgf, iset), la_min=la_min(iset), &
923 lb_max=0, zetb=0.0_dp, lb_min=0, &
924 ra=ra, rab=[0.0_dp, 0.0_dp, 0.0_dp], &
925 rsgrid=rs_grid, &
926 hab=hab, pab=pab, o1=na1, o2=0, radius=radius, &
927 calculate_forces=calculate_forces, &
928 force_a=force_a, force_b=force_b, &
929 use_virial=use_virial, &
930 my_virial_a=my_virial_a, my_virial_b=my_virial_b)
931 END IF
932 END IF
933 END DO
934 ! contract hab
935 CALL dgemm("T", "N", nsgf_seta(iset), 1, ncoa, 1.0_dp, sphi_a(1, sgfa), &
936 SIZE(sphi_a, 1), hab(1, 1), SIZE(hab, 1), 0.0_dp, work_i(1, 1), SIZE(work_i, 1))
937
938 int_res(ikind)%v_int(iatom, sgfa:sgfa - 1 + nsgf_seta(iset)) = &
939 int_res(ikind)%v_int(iatom, sgfa:sgfa - 1 + nsgf_seta(iset)) + work_i(1:nsgf_seta(iset), 1)
940 END IF
941 END DO
942 !
943 IF (calculate_forces) THEN
944 int_res(ikind)%v_dfdr(iatom, :) = int_res(ikind)%v_dfdr(iatom, :) + force_a(:)
945 IF (use_virial) THEN
946 virial%pv_lrigpw = virial%pv_lrigpw + my_virial_a
947 virial%pv_virial = virial%pv_virial + my_virial_a
948 END IF
949 END IF
950
951 END DO
952
953 IF (calculate_forces) DEALLOCATE (work_f)
954 DEALLOCATE (work_i, map_it)
955 DEALLOCATE (hab, pab)
956 END DO
957
958 CALL timestop(handle)
959
960 END SUBROUTINE integrate_v_rspace_one_center
961
962! **************************************************************************************************
963!> \brief computes integrals of product of v_rspace times the diagonal block basis functions
964!> required for LRIGPW with exact 1c terms
965!> \param v_rspace ...
966!> \param ksmat ...
967!> \param pmat ...
968!> \param qs_env ...
969!> \param calculate_forces ...
970!> \param basis_type ...
971!> \author JGH
972! **************************************************************************************************
973 SUBROUTINE integrate_v_rspace_diagonal(v_rspace, ksmat, pmat, qs_env, calculate_forces, basis_type)
974 TYPE(pw_r3d_rs_type), INTENT(IN) :: v_rspace
975 TYPE(dbcsr_type), INTENT(INOUT) :: ksmat, pmat
976 TYPE(qs_environment_type), POINTER :: qs_env
977 LOGICAL, INTENT(IN) :: calculate_forces
978 CHARACTER(len=*), INTENT(IN) :: basis_type
979
980 CHARACTER(len=*), PARAMETER :: routinen = 'integrate_v_rspace_diagonal'
981
982 INTEGER :: atom_a, group_size, handle, iatom, igrid_level, ikind, ipgf, iset, jpgf, jset, &
983 m1, maxco, maxsgf_set, my_pos, na1, na2, natom_of_kind, nb1, nb2, ncoa, ncob, nkind, &
984 nseta, nsgfa, offset, sgfa, sgfb
985 INTEGER, DIMENSION(:), POINTER :: atom_list, la_max, la_min, npgfa, &
986 nsgf_seta
987 INTEGER, DIMENSION(:, :), POINTER :: first_sgfa
988 LOGICAL :: found, use_virial
989 LOGICAL, ALLOCATABLE, DIMENSION(:, :) :: map_it2
990 REAL(kind=dp) :: eps_rho_rspace, radius, zetp
991 REAL(kind=dp), DIMENSION(3) :: force_a, force_b, ra
992 REAL(kind=dp), DIMENSION(3, 3) :: my_virial_a, my_virial_b
993 REAL(kind=dp), DIMENSION(:), POINTER :: set_radius_a
994 REAL(kind=dp), DIMENSION(:, :), POINTER :: h_block, hab, hmat, p_block, pab, pblk, &
995 rpgfa, sphi_a, work, zeta
996 TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
997 TYPE(cell_type), POINTER :: cell
998 TYPE(dft_control_type), POINTER :: dft_control
999 TYPE(gridlevel_info_type), POINTER :: gridlevel_info
1000 TYPE(gto_basis_set_type), POINTER :: orb_basis_set
1001 TYPE(mp_para_env_type), POINTER :: para_env
1002 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
1003 TYPE(pw_env_type), POINTER :: pw_env
1004 TYPE(qs_force_type), DIMENSION(:), POINTER :: force
1005 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
1006 TYPE(realspace_grid_type), DIMENSION(:), POINTER :: rs_v
1007 TYPE(realspace_grid_type), POINTER :: rs_grid
1008 TYPE(virial_type), POINTER :: virial
1009
1010 CALL timeset(routinen, handle)
1011
1012 CALL get_qs_env(qs_env=qs_env, pw_env=pw_env)
1013 CALL pw_env_get(pw_env, rs_grids=rs_v)
1014 CALL potential_pw2rs(rs_v, v_rspace, pw_env)
1015
1016 gridlevel_info => pw_env%gridlevel_info
1017
1018 CALL get_qs_env(qs_env=qs_env, &
1019 atomic_kind_set=atomic_kind_set, &
1020 qs_kind_set=qs_kind_set, &
1021 cell=cell, &
1022 dft_control=dft_control, &
1023 nkind=nkind, &
1024 particle_set=particle_set, &
1025 force=force, &
1026 virial=virial, &
1027 para_env=para_env)
1028
1029 eps_rho_rspace = dft_control%qs_control%eps_rho_rspace
1030 use_virial = virial%pv_availability .AND. (.NOT. virial%pv_numer)
1031
1032 offset = 0
1033 my_pos = v_rspace%pw_grid%para%group%mepos
1034 group_size = v_rspace%pw_grid%para%group%num_pe
1035
1036 DO ikind = 1, nkind
1037
1038 CALL get_atomic_kind(atomic_kind_set(ikind), natom=natom_of_kind, atom_list=atom_list)
1039 CALL get_qs_kind(qs_kind_set(ikind), basis_set=orb_basis_set, basis_type=basis_type)
1040 CALL get_gto_basis_set(gto_basis_set=orb_basis_set, &
1041 lmax=la_max, lmin=la_min, maxco=maxco, maxsgf_set=maxsgf_set, &
1042 npgf=npgfa, nset=nseta, nsgf_set=nsgf_seta, nsgf=nsgfa, &
1043 first_sgf=first_sgfa, pgf_radius=rpgfa, set_radius=set_radius_a, &
1044 sphi=sphi_a, zet=zeta)
1045
1046 ALLOCATE (hab(maxco, maxco), work(maxco, maxsgf_set), hmat(nsgfa, nsgfa))
1047 IF (calculate_forces) ALLOCATE (pab(maxco, maxco), pblk(nsgfa, nsgfa))
1048
1049 DO iatom = 1, natom_of_kind
1050 atom_a = atom_list(iatom)
1051 ra(:) = pbc(particle_set(atom_a)%r, cell)
1052 hmat = 0.0_dp
1053 IF (calculate_forces) THEN
1054 CALL dbcsr_get_block_p(matrix=pmat, row=atom_a, col=atom_a, block=p_block, found=found)
1055 IF (found) THEN
1056 pblk(1:nsgfa, 1:nsgfa) = p_block(1:nsgfa, 1:nsgfa)
1057 ELSE
1058 pblk = 0.0_dp
1059 END IF
1060 CALL para_env%sum(pblk)
1061 force_a(:) = 0._dp
1062 force_b(:) = 0._dp
1063 IF (use_virial) THEN
1064 my_virial_a = 0.0_dp
1065 my_virial_b = 0.0_dp
1066 END IF
1067 END IF
1068 m1 = maxval(npgfa(1:nseta))
1069 ALLOCATE (map_it2(m1, m1))
1070 DO iset = 1, nseta
1071 sgfa = first_sgfa(1, iset)
1072 ncoa = npgfa(iset)*ncoset(la_max(iset))
1073 DO jset = 1, nseta
1074 sgfb = first_sgfa(1, jset)
1075 ncob = npgfa(jset)*ncoset(la_max(jset))
1076 !
1077 map_it2 = .false.
1078 DO ipgf = 1, npgfa(iset)
1079 DO jpgf = 1, npgfa(jset)
1080 zetp = zeta(ipgf, iset) + zeta(jpgf, jset)
1081 igrid_level = gaussian_gridlevel(gridlevel_info, zetp)
1082 rs_grid => rs_v(igrid_level)
1083 map_it2(ipgf, jpgf) = map_gaussian_here(rs_grid, cell%h_inv, ra, offset, group_size, my_pos)
1084 END DO
1085 END DO
1086 offset = offset + 1
1087 !
1088 IF (any(map_it2(1:npgfa(iset), 1:npgfa(jset)))) THEN
1089 hab(:, :) = 0._dp
1090 IF (calculate_forces) THEN
1091 CALL dgemm("N", "N", ncoa, nsgf_seta(jset), nsgf_seta(iset), &
1092 1.0_dp, sphi_a(1, sgfa), SIZE(sphi_a, 1), &
1093 pblk(sgfa, sgfb), SIZE(pblk, 1), &
1094 0.0_dp, work(1, 1), SIZE(work, 1))
1095 CALL dgemm("N", "T", ncoa, ncob, nsgf_seta(jset), &
1096 1.0_dp, work(1, 1), SIZE(work, 1), &
1097 sphi_a(1, sgfb), SIZE(sphi_a, 1), &
1098 0.0_dp, pab(1, 1), SIZE(pab, 1))
1099 END IF
1100
1101 DO ipgf = 1, npgfa(iset)
1102 na1 = (ipgf - 1)*ncoset(la_max(iset))
1103 na2 = ipgf*ncoset(la_max(iset))
1104 DO jpgf = 1, npgfa(jset)
1105 nb1 = (jpgf - 1)*ncoset(la_max(jset))
1106 nb2 = jpgf*ncoset(la_max(jset))
1107 zetp = zeta(ipgf, iset) + zeta(jpgf, jset)
1108 igrid_level = gaussian_gridlevel(gridlevel_info, zetp)
1109 rs_grid => rs_v(igrid_level)
1110
1111 radius = exp_radius_very_extended(la_min=la_min(iset), la_max=la_max(iset), &
1112 lb_min=la_min(jset), lb_max=la_max(jset), &
1113 ra=ra, rb=ra, rp=ra, &
1114 zetp=zetp, eps=eps_rho_rspace, &
1115 prefactor=1.0_dp, cutoff=1.0_dp)
1116
1117 IF (map_it2(ipgf, jpgf)) THEN
1118 IF (calculate_forces) THEN
1119 CALL integrate_pgf_product( &
1120 la_max(iset), zeta(ipgf, iset), la_min(iset), &
1121 la_max(jset), zeta(jpgf, jset), la_min(jset), &
1122 ra=ra, rab=[0.0_dp, 0.0_dp, 0.0_dp], &
1123 rsgrid=rs_v(igrid_level), &
1124 hab=hab, pab=pab, o1=na1, o2=nb1, &
1125 radius=radius, &
1126 calculate_forces=.true., &
1127 force_a=force_a, force_b=force_b, &
1128 use_virial=use_virial, my_virial_a=my_virial_a, my_virial_b=my_virial_b)
1129 ELSE
1130 CALL integrate_pgf_product( &
1131 la_max(iset), zeta(ipgf, iset), la_min(iset), &
1132 la_max(jset), zeta(jpgf, jset), la_min(jset), &
1133 ra=ra, rab=[0.0_dp, 0.0_dp, 0.0_dp], &
1134 rsgrid=rs_v(igrid_level), &
1135 hab=hab, o1=na1, o2=nb1, &
1136 radius=radius, &
1137 calculate_forces=.false.)
1138 END IF
1139 END IF
1140 END DO
1141 END DO
1142 ! contract hab
1143 CALL dgemm("N", "N", ncoa, nsgf_seta(jset), ncob, &
1144 1.0_dp, hab(1, 1), SIZE(hab, 1), &
1145 sphi_a(1, sgfb), SIZE(sphi_a, 1), &
1146 0.0_dp, work(1, 1), SIZE(work, 1))
1147 CALL dgemm("T", "N", nsgf_seta(iset), nsgf_seta(jset), ncoa, &
1148 1.0_dp, sphi_a(1, sgfa), SIZE(sphi_a, 1), &
1149 work(1, 1), SIZE(work, 1), &
1150 1.0_dp, hmat(sgfa, sgfb), SIZE(hmat, 1))
1151 END IF
1152 END DO
1153 END DO
1154 DEALLOCATE (map_it2)
1155 ! update KS matrix
1156 CALL para_env%sum(hmat)
1157 CALL dbcsr_get_block_p(matrix=ksmat, row=atom_a, col=atom_a, block=h_block, found=found)
1158 IF (found) THEN
1159 h_block(1:nsgfa, 1:nsgfa) = h_block(1:nsgfa, 1:nsgfa) + hmat(1:nsgfa, 1:nsgfa)
1160 END IF
1161 IF (calculate_forces) THEN
1162 force(ikind)%rho_elec(:, iatom) = force(ikind)%rho_elec(:, iatom) + 2.0_dp*force_a(:)
1163 IF (use_virial) THEN
1164 IF (use_virial .AND. calculate_forces) THEN
1165 virial%pv_lrigpw = virial%pv_lrigpw + 2.0_dp*my_virial_a
1166 virial%pv_virial = virial%pv_virial + 2.0_dp*my_virial_a
1167 END IF
1168 END IF
1169 END IF
1170 END DO
1171 DEALLOCATE (hab, work, hmat)
1172 IF (calculate_forces) DEALLOCATE (pab, pblk)
1173 END DO
1174
1175 CALL timestop(handle)
1176
1177 END SUBROUTINE integrate_v_rspace_diagonal
1178
1179! **************************************************************************************************
1180!> \brief computes integrals of product of v_rspace times a basis function (vector function)
1181!> and possible forces
1182!> \param qs_env ...
1183!> \param v_rspace ...
1184!> \param f_coef ...
1185!> \param f_integral ...
1186!> \param calculate_forces ...
1187!> \param basis_type ...
1188!> \author JGH [8.2024]
1189! **************************************************************************************************
1190 SUBROUTINE integrate_function(qs_env, v_rspace, f_coef, f_integral, calculate_forces, basis_type)
1191 TYPE(qs_environment_type), POINTER :: qs_env
1192 TYPE(pw_r3d_rs_type), INTENT(IN) :: v_rspace
1193 REAL(kind=dp), DIMENSION(:), INTENT(IN) :: f_coef
1194 REAL(kind=dp), DIMENSION(:), INTENT(OUT) :: f_integral
1195 LOGICAL, INTENT(IN) :: calculate_forces
1196 CHARACTER(len=*), INTENT(IN) :: basis_type
1197
1198 CHARACTER(len=*), PARAMETER :: routinen = 'integrate_function'
1199
1200 INTEGER :: atom_a, group_size, handle, i, iatom, igrid_level, ikind, ipgf, iset, maxco, &
1201 maxsgf_set, my_pos, na1, natom, ncoa, nkind, nseta, offset, sgfa
1202 INTEGER, ALLOCATABLE, DIMENSION(:) :: atom_of_kind
1203 INTEGER, DIMENSION(:), POINTER :: la_max, la_min, npgfa, nsgfa
1204 INTEGER, DIMENSION(:, :), POINTER :: first_sgfa
1205 LOGICAL :: use_virial
1206 REAL(kind=dp) :: eps_rho_rspace, radius
1207 REAL(kind=dp), DIMENSION(3) :: force_a, force_b, ra
1208 REAL(kind=dp), DIMENSION(3, 3) :: my_virial_a, my_virial_b
1209 REAL(kind=dp), DIMENSION(:, :), POINTER :: hab, pab, sphi_a, work, zeta
1210 TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
1211 TYPE(cell_type), POINTER :: cell
1212 TYPE(dft_control_type), POINTER :: dft_control
1213 TYPE(gridlevel_info_type), POINTER :: gridlevel_info
1214 TYPE(gto_basis_set_type), POINTER :: orb_basis_set
1215 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
1216 TYPE(pw_env_type), POINTER :: pw_env
1217 TYPE(qs_force_type), DIMENSION(:), POINTER :: force
1218 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
1219 TYPE(realspace_grid_type), DIMENSION(:), POINTER :: rs_v
1220 TYPE(realspace_grid_type), POINTER :: rs_grid
1221 TYPE(virial_type), POINTER :: virial
1222
1223 CALL timeset(routinen, handle)
1224
1225 CALL get_qs_env(qs_env=qs_env, pw_env=pw_env)
1226 gridlevel_info => pw_env%gridlevel_info
1227
1228 CALL pw_env_get(pw_env, rs_grids=rs_v)
1229 CALL potential_pw2rs(rs_v, v_rspace, pw_env)
1230
1231 CALL get_qs_env(qs_env=qs_env, &
1232 atomic_kind_set=atomic_kind_set, &
1233 qs_kind_set=qs_kind_set, &
1234 force=force, &
1235 cell=cell, &
1236 dft_control=dft_control, &
1237 nkind=nkind, &
1238 natom=natom, &
1239 particle_set=particle_set, &
1240 pw_env=pw_env, &
1241 virial=virial)
1242 CALL get_atomic_kind_set(atomic_kind_set=atomic_kind_set, atom_of_kind=atom_of_kind)
1243
1244 use_virial = virial%pv_availability .AND. (.NOT. virial%pv_numer)
1245 IF (use_virial) THEN
1246 cpabort("Virial NYA")
1247 END IF
1248
1249 eps_rho_rspace = dft_control%qs_control%eps_rho_rspace
1250
1251 CALL get_qs_kind_set(qs_kind_set, &
1252 maxco=maxco, maxsgf_set=maxsgf_set, basis_type=basis_type)
1253 ALLOCATE (hab(maxco, 1), pab(maxco, 1), work(maxco, 1))
1254
1255 offset = 0
1256 my_pos = v_rspace%pw_grid%para%group%mepos
1257 group_size = v_rspace%pw_grid%para%group%num_pe
1258
1259 DO iatom = 1, natom
1260 ikind = particle_set(iatom)%atomic_kind%kind_number
1261 atom_a = atom_of_kind(iatom)
1262 CALL get_qs_kind(qs_kind_set(ikind), basis_set=orb_basis_set, basis_type=basis_type)
1263 CALL get_gto_basis_set(gto_basis_set=orb_basis_set, &
1264 first_sgf=first_sgfa, &
1265 lmax=la_max, &
1266 lmin=la_min, &
1267 npgf=npgfa, &
1268 nset=nseta, &
1269 nsgf_set=nsgfa, &
1270 sphi=sphi_a, &
1271 zet=zeta)
1272 ra(:) = pbc(particle_set(iatom)%r, cell)
1273
1274 force_a(:) = 0._dp
1275 force_b(:) = 0._dp
1276 my_virial_a(:, :) = 0._dp
1277 my_virial_b(:, :) = 0._dp
1278
1279 DO iset = 1, nseta
1280
1281 ncoa = npgfa(iset)*ncoset(la_max(iset))
1282 sgfa = first_sgfa(1, iset)
1283
1284 hab = 0._dp
1285 pab = 0._dp
1286
1287 DO i = 1, nsgfa(iset)
1288 work(i, 1) = f_coef(offset + i)
1289 END DO
1290
1291 CALL dgemm("N", "N", ncoa, 1, nsgfa(iset), &
1292 1.0_dp, sphi_a(1, sgfa), SIZE(sphi_a, 1), &
1293 work(1, 1), SIZE(work, 1), &
1294 0.0_dp, pab(1, 1), SIZE(pab, 1))
1295
1296 DO ipgf = 1, npgfa(iset)
1297
1298 na1 = (ipgf - 1)*ncoset(la_max(iset))
1299
1300 igrid_level = gaussian_gridlevel(gridlevel_info, zeta(ipgf, iset))
1301 rs_grid => rs_v(igrid_level)
1302
1303 IF (map_gaussian_here(rs_grid, cell%h_inv, ra, offset, group_size, my_pos)) THEN
1304 radius = exp_radius_very_extended(la_min=la_min(iset), la_max=la_max(iset), &
1305 lb_min=0, lb_max=0, ra=ra, rb=ra, rp=ra, &
1306 zetp=zeta(ipgf, iset), eps=eps_rho_rspace, &
1307 prefactor=1.0_dp, cutoff=1.0_dp)
1308
1309 IF (calculate_forces) THEN
1310 CALL integrate_pgf_product(la_max=la_max(iset), &
1311 zeta=zeta(ipgf, iset), la_min=la_min(iset), &
1312 lb_max=0, zetb=0.0_dp, lb_min=0, &
1313 ra=ra, rab=[0.0_dp, 0.0_dp, 0.0_dp], &
1314 rsgrid=rs_grid, &
1315 hab=hab, pab=pab, o1=na1, o2=0, radius=radius, &
1316 calculate_forces=.true., &
1317 force_a=force_a, force_b=force_b, &
1318 use_virial=use_virial, &
1319 my_virial_a=my_virial_a, my_virial_b=my_virial_b)
1320 ELSE
1321 CALL integrate_pgf_product(la_max=la_max(iset), &
1322 zeta=zeta(ipgf, iset), la_min=la_min(iset), &
1323 lb_max=0, zetb=0.0_dp, lb_min=0, &
1324 ra=ra, rab=[0.0_dp, 0.0_dp, 0.0_dp], &
1325 rsgrid=rs_grid, &
1326 hab=hab, o1=na1, o2=0, radius=radius, &
1327 calculate_forces=.false.)
1328 END IF
1329
1330 END IF
1331
1332 END DO
1333 !
1334 CALL dgemm("T", "N", nsgfa(iset), 1, ncoa, 1.0_dp, sphi_a(1, sgfa), &
1335 SIZE(sphi_a, 1), hab(1, 1), SIZE(hab, 1), 0.0_dp, work(1, 1), SIZE(work, 1))
1336 DO i = 1, nsgfa(iset)
1337 f_integral(offset + i) = work(i, 1)
1338 END DO
1339
1340 offset = offset + nsgfa(iset)
1341
1342 END DO
1343
1344 IF (calculate_forces) THEN
1345 force(ikind)%rho_elec(:, atom_a) = force(ikind)%rho_elec(:, atom_a) + force_a(:)
1346 IF (use_virial) THEN
1347 virial%pv_virial = virial%pv_virial + my_virial_a
1348 END IF
1349 END IF
1350
1351 END DO
1352
1353 DEALLOCATE (hab, pab, work)
1354
1355 CALL timestop(handle)
1356
1357 END SUBROUTINE integrate_function
1358
1359END 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:218
cell_types::pbc
Definition cell_types.F:97
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)
...
Definition basis_set_types.F:631
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:692
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:277
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:530
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:1191
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:649
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:974
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:486
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:767
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:286
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:694
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:2144
realspace_grid_types::rs_grid_zero
subroutine, public rs_grid_zero(rs)
Initialize grid to zero.
Definition realspace_grid_types.F:1905
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:633
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:743
particle_types::particle_type
Definition particle_types.F:35
pw_env_types::pw_env_type
contained for different pw related things
Definition pw_env_types.F:70
pw_types::pw_r3d_rs_type
Definition pw_types.F:62
qs_environment_types::qs_environment_type
Definition qs_environment_types.F:220
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