(git:b17b328)
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-2025 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! **************************************************************************************************
39 USE cell_types, ONLY: cell_type,&
40 pbc
49 USE kinds, ONLY: dp
53 USE orbital_pointers, ONLY: coset,&
54 ncoset
56 USE pw_env_types, ONLY: pw_env_get,&
58 USE pw_types, ONLY: pw_r3d_rs_type
62 USE qs_kind_types, ONLY: get_qs_kind,&
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
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, maxco, &
777 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
860 DO iatom = 1, natom_of_kind
861
862 atom_a = atom_list(iatom)
863 IF (PRESENT(atomlist)) THEN
864 IF (atomlist(atom_a) == 0) cycle
865 END IF
866 ra(:) = pbc(particle_set(atom_a)%r, cell)
867 force_a(:) = 0._dp
868 force_b(:) = 0._dp
869 my_virial_a(:, :) = 0._dp
870 my_virial_b(:, :) = 0._dp
871
872 m1 = maxval(npgfa(1:nseta))
873 ALLOCATE (map_it(m1))
874
875 DO iset = 1, nseta
876 !
877 map_it = .false.
878 DO ipgf = 1, npgfa(iset)
879 igrid_level = gaussian_gridlevel(gridlevel_info, zeta(ipgf, iset))
880 rs_grid => rs_v(igrid_level)
881 map_it(ipgf) = map_gaussian_here(rs_grid, cell%h_inv, ra, offset, group_size, my_pos)
882 END DO
883 offset = offset + 1
884 !
885 IF (any(map_it(1:npgfa(iset)))) THEN
886 sgfa = first_sgfa(1, iset)
887 ncoa = npgfa(iset)*ncoset(la_max(iset))
888 hab(:, 1) = 0._dp
889 ALLOCATE (work_i(nsgf_seta(iset), 1))
890 work_i = 0.0_dp
891
892 ! get fit coefficients for forces
893 IF (calculate_forces) THEN
894 m1 = sgfa + nsgf_seta(iset) - 1
895 ALLOCATE (work_f(nsgf_seta(iset), 1))
896 work_f(1:nsgf_seta(iset), 1) = int_res(ikind)%acoef(iatom, sgfa:m1)
897 CALL dgemm("N", "N", ncoa, 1, nsgf_seta(iset), 1.0_dp, sphi_a(1, sgfa), &
898 SIZE(sphi_a, 1), work_f(1, 1), SIZE(work_f, 1), 0.0_dp, pab(1, 1), &
899 SIZE(pab, 1))
900 DEALLOCATE (work_f)
901 END IF
902
903 DO ipgf = 1, npgfa(iset)
904 na1 = (ipgf - 1)*ncoset(la_max(iset))
905 igrid_level = gaussian_gridlevel(gridlevel_info, zeta(ipgf, iset))
906 rs_grid => rs_v(igrid_level)
907
908 radius = exp_radius_very_extended(la_min=la_min(iset), la_max=la_max(iset), &
909 lb_min=0, lb_max=0, ra=ra, rb=ra, rp=ra, &
910 zetp=zeta(ipgf, iset), eps=eps_rho_rspace, &
911 prefactor=1.0_dp, cutoff=1.0_dp)
912
913 IF (map_it(ipgf)) THEN
914 IF (.NOT. calculate_forces) THEN
915 CALL integrate_pgf_product(la_max=la_max(iset), &
916 zeta=zeta(ipgf, iset), la_min=la_min(iset), &
917 lb_max=0, zetb=0.0_dp, lb_min=0, &
918 ra=ra, rab=(/0.0_dp, 0.0_dp, 0.0_dp/), &
919 rsgrid=rs_grid, &
920 hab=hab, o1=na1, o2=0, radius=radius, &
921 calculate_forces=calculate_forces)
922 ELSE
923 CALL integrate_pgf_product(la_max=la_max(iset), &
924 zeta=zeta(ipgf, iset), la_min=la_min(iset), &
925 lb_max=0, zetb=0.0_dp, lb_min=0, &
926 ra=ra, rab=(/0.0_dp, 0.0_dp, 0.0_dp/), &
927 rsgrid=rs_grid, &
928 hab=hab, pab=pab, o1=na1, o2=0, radius=radius, &
929 calculate_forces=calculate_forces, &
930 force_a=force_a, force_b=force_b, &
931 use_virial=use_virial, &
932 my_virial_a=my_virial_a, my_virial_b=my_virial_b)
933 END IF
934 END IF
935 END DO
936 ! contract hab
937 CALL dgemm("T", "N", nsgf_seta(iset), 1, ncoa, 1.0_dp, sphi_a(1, sgfa), &
938 SIZE(sphi_a, 1), hab(1, 1), SIZE(hab, 1), 0.0_dp, work_i(1, 1), SIZE(work_i, 1))
939
940 int_res(ikind)%v_int(iatom, sgfa:sgfa - 1 + nsgf_seta(iset)) = &
941 int_res(ikind)%v_int(iatom, sgfa:sgfa - 1 + nsgf_seta(iset)) + work_i(1:nsgf_seta(iset), 1)
942 DEALLOCATE (work_i)
943 END IF
944 END DO
945 !
946 IF (calculate_forces) THEN
947 int_res(ikind)%v_dfdr(iatom, :) = int_res(ikind)%v_dfdr(iatom, :) + force_a(:)
948 IF (use_virial) THEN
949 virial%pv_lrigpw = virial%pv_lrigpw + my_virial_a
950 virial%pv_virial = virial%pv_virial + my_virial_a
951 END IF
952 END IF
953
954 DEALLOCATE (map_it)
955
956 END DO
957
958 DEALLOCATE (hab, pab)
959 END DO
960
961 CALL timestop(handle)
962
963 END SUBROUTINE integrate_v_rspace_one_center
964
965! **************************************************************************************************
966!> \brief computes integrals of product of v_rspace times the diagonal block basis functions
967!> required for LRIGPW with exact 1c terms
968!> \param v_rspace ...
969!> \param ksmat ...
970!> \param pmat ...
971!> \param qs_env ...
972!> \param calculate_forces ...
973!> \param basis_type ...
974!> \author JGH
975! **************************************************************************************************
976 SUBROUTINE integrate_v_rspace_diagonal(v_rspace, ksmat, pmat, qs_env, calculate_forces, basis_type)
977 TYPE(pw_r3d_rs_type), INTENT(IN) :: v_rspace
978 TYPE(dbcsr_type), INTENT(INOUT) :: ksmat, pmat
979 TYPE(qs_environment_type), POINTER :: qs_env
980 LOGICAL, INTENT(IN) :: calculate_forces
981 CHARACTER(len=*), INTENT(IN) :: basis_type
982
983 CHARACTER(len=*), PARAMETER :: routinen = 'integrate_v_rspace_diagonal'
984
985 INTEGER :: atom_a, group_size, handle, iatom, igrid_level, ikind, ipgf, iset, jpgf, jset, &
986 m1, maxco, maxsgf_set, my_pos, na1, na2, natom_of_kind, nb1, nb2, ncoa, ncob, nkind, &
987 nseta, nsgfa, offset, sgfa, sgfb
988 INTEGER, DIMENSION(:), POINTER :: atom_list, la_max, la_min, npgfa, &
989 nsgf_seta
990 INTEGER, DIMENSION(:, :), POINTER :: first_sgfa
991 LOGICAL :: found, use_virial
992 LOGICAL, ALLOCATABLE, DIMENSION(:, :) :: map_it2
993 REAL(kind=dp) :: eps_rho_rspace, radius, zetp
994 REAL(kind=dp), DIMENSION(3) :: force_a, force_b, ra
995 REAL(kind=dp), DIMENSION(3, 3) :: my_virial_a, my_virial_b
996 REAL(kind=dp), DIMENSION(:), POINTER :: set_radius_a
997 REAL(kind=dp), DIMENSION(:, :), POINTER :: h_block, hab, hmat, p_block, pab, pblk, &
998 rpgfa, sphi_a, work, zeta
999 TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
1000 TYPE(cell_type), POINTER :: cell
1001 TYPE(dft_control_type), POINTER :: dft_control
1002 TYPE(gridlevel_info_type), POINTER :: gridlevel_info
1003 TYPE(gto_basis_set_type), POINTER :: orb_basis_set
1004 TYPE(mp_para_env_type), POINTER :: para_env
1005 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
1006 TYPE(pw_env_type), POINTER :: pw_env
1007 TYPE(qs_force_type), DIMENSION(:), POINTER :: force
1008 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
1009 TYPE(realspace_grid_type), DIMENSION(:), POINTER :: rs_v
1010 TYPE(realspace_grid_type), POINTER :: rs_grid
1011 TYPE(virial_type), POINTER :: virial
1012
1013 CALL timeset(routinen, handle)
1014
1015 CALL get_qs_env(qs_env=qs_env, pw_env=pw_env)
1016 CALL pw_env_get(pw_env, rs_grids=rs_v)
1017 CALL potential_pw2rs(rs_v, v_rspace, pw_env)
1018
1019 gridlevel_info => pw_env%gridlevel_info
1020
1021 CALL get_qs_env(qs_env=qs_env, &
1022 atomic_kind_set=atomic_kind_set, &
1023 qs_kind_set=qs_kind_set, &
1024 cell=cell, &
1025 dft_control=dft_control, &
1026 nkind=nkind, &
1027 particle_set=particle_set, &
1028 force=force, &
1029 virial=virial, &
1030 para_env=para_env)
1031
1032 eps_rho_rspace = dft_control%qs_control%eps_rho_rspace
1033 use_virial = virial%pv_availability .AND. (.NOT. virial%pv_numer)
1034
1035 offset = 0
1036 my_pos = v_rspace%pw_grid%para%group%mepos
1037 group_size = v_rspace%pw_grid%para%group%num_pe
1038
1039 DO ikind = 1, nkind
1040
1041 CALL get_atomic_kind(atomic_kind_set(ikind), natom=natom_of_kind, atom_list=atom_list)
1042 CALL get_qs_kind(qs_kind_set(ikind), basis_set=orb_basis_set, basis_type=basis_type)
1043 CALL get_gto_basis_set(gto_basis_set=orb_basis_set, &
1044 lmax=la_max, lmin=la_min, maxco=maxco, maxsgf_set=maxsgf_set, &
1045 npgf=npgfa, nset=nseta, nsgf_set=nsgf_seta, nsgf=nsgfa, &
1046 first_sgf=first_sgfa, pgf_radius=rpgfa, set_radius=set_radius_a, &
1047 sphi=sphi_a, zet=zeta)
1048
1049 ALLOCATE (hab(maxco, maxco), work(maxco, maxsgf_set), hmat(nsgfa, nsgfa))
1050 IF (calculate_forces) ALLOCATE (pab(maxco, maxco), pblk(nsgfa, nsgfa))
1051
1052 DO iatom = 1, natom_of_kind
1053 atom_a = atom_list(iatom)
1054 ra(:) = pbc(particle_set(atom_a)%r, cell)
1055 hmat = 0.0_dp
1056 IF (calculate_forces) THEN
1057 CALL dbcsr_get_block_p(matrix=pmat, row=atom_a, col=atom_a, block=p_block, found=found)
1058 IF (found) THEN
1059 pblk(1:nsgfa, 1:nsgfa) = p_block(1:nsgfa, 1:nsgfa)
1060 ELSE
1061 pblk = 0.0_dp
1062 END IF
1063 CALL para_env%sum(pblk)
1064 force_a(:) = 0._dp
1065 force_b(:) = 0._dp
1066 IF (use_virial) THEN
1067 my_virial_a = 0.0_dp
1068 my_virial_b = 0.0_dp
1069 END IF
1070 END IF
1071 m1 = maxval(npgfa(1:nseta))
1072 ALLOCATE (map_it2(m1, m1))
1073 DO iset = 1, nseta
1074 sgfa = first_sgfa(1, iset)
1075 ncoa = npgfa(iset)*ncoset(la_max(iset))
1076 DO jset = 1, nseta
1077 sgfb = first_sgfa(1, jset)
1078 ncob = npgfa(jset)*ncoset(la_max(jset))
1079 !
1080 map_it2 = .false.
1081 DO ipgf = 1, npgfa(iset)
1082 DO jpgf = 1, npgfa(jset)
1083 zetp = zeta(ipgf, iset) + zeta(jpgf, jset)
1084 igrid_level = gaussian_gridlevel(gridlevel_info, zetp)
1085 rs_grid => rs_v(igrid_level)
1086 map_it2(ipgf, jpgf) = map_gaussian_here(rs_grid, cell%h_inv, ra, offset, group_size, my_pos)
1087 END DO
1088 END DO
1089 offset = offset + 1
1090 !
1091 IF (any(map_it2(1:npgfa(iset), 1:npgfa(jset)))) THEN
1092 hab(:, :) = 0._dp
1093 IF (calculate_forces) THEN
1094 CALL dgemm("N", "N", ncoa, nsgf_seta(jset), nsgf_seta(iset), &
1095 1.0_dp, sphi_a(1, sgfa), SIZE(sphi_a, 1), &
1096 pblk(sgfa, sgfb), SIZE(pblk, 1), &
1097 0.0_dp, work(1, 1), SIZE(work, 1))
1098 CALL dgemm("N", "T", ncoa, ncob, nsgf_seta(jset), &
1099 1.0_dp, work(1, 1), SIZE(work, 1), &
1100 sphi_a(1, sgfb), SIZE(sphi_a, 1), &
1101 0.0_dp, pab(1, 1), SIZE(pab, 1))
1102 END IF
1103
1104 DO ipgf = 1, npgfa(iset)
1105 na1 = (ipgf - 1)*ncoset(la_max(iset))
1106 na2 = ipgf*ncoset(la_max(iset))
1107 DO jpgf = 1, npgfa(jset)
1108 nb1 = (jpgf - 1)*ncoset(la_max(jset))
1109 nb2 = jpgf*ncoset(la_max(jset))
1110 zetp = zeta(ipgf, iset) + zeta(jpgf, jset)
1111 igrid_level = gaussian_gridlevel(gridlevel_info, zetp)
1112 rs_grid => rs_v(igrid_level)
1113
1114 radius = exp_radius_very_extended(la_min=la_min(iset), la_max=la_max(iset), &
1115 lb_min=la_min(jset), lb_max=la_max(jset), &
1116 ra=ra, rb=ra, rp=ra, &
1117 zetp=zetp, eps=eps_rho_rspace, &
1118 prefactor=1.0_dp, cutoff=1.0_dp)
1119
1120 IF (map_it2(ipgf, jpgf)) THEN
1121 IF (calculate_forces) THEN
1122 CALL integrate_pgf_product( &
1123 la_max(iset), zeta(ipgf, iset), la_min(iset), &
1124 la_max(jset), zeta(jpgf, jset), la_min(jset), &
1125 ra=ra, rab=(/0.0_dp, 0.0_dp, 0.0_dp/), &
1126 rsgrid=rs_v(igrid_level), &
1127 hab=hab, pab=pab, o1=na1, o2=nb1, &
1128 radius=radius, &
1129 calculate_forces=.true., &
1130 force_a=force_a, force_b=force_b, &
1131 use_virial=use_virial, my_virial_a=my_virial_a, my_virial_b=my_virial_b)
1132 ELSE
1133 CALL integrate_pgf_product( &
1134 la_max(iset), zeta(ipgf, iset), la_min(iset), &
1135 la_max(jset), zeta(jpgf, jset), la_min(jset), &
1136 ra=ra, rab=(/0.0_dp, 0.0_dp, 0.0_dp/), &
1137 rsgrid=rs_v(igrid_level), &
1138 hab=hab, o1=na1, o2=nb1, &
1139 radius=radius, &
1140 calculate_forces=.false.)
1141 END IF
1142 END IF
1143 END DO
1144 END DO
1145 ! contract hab
1146 CALL dgemm("N", "N", ncoa, nsgf_seta(jset), ncob, &
1147 1.0_dp, hab(1, 1), SIZE(hab, 1), &
1148 sphi_a(1, sgfb), SIZE(sphi_a, 1), &
1149 0.0_dp, work(1, 1), SIZE(work, 1))
1150 CALL dgemm("T", "N", nsgf_seta(iset), nsgf_seta(jset), ncoa, &
1151 1.0_dp, sphi_a(1, sgfa), SIZE(sphi_a, 1), &
1152 work(1, 1), SIZE(work, 1), &
1153 1.0_dp, hmat(sgfa, sgfb), SIZE(hmat, 1))
1154 END IF
1155 END DO
1156 END DO
1157 DEALLOCATE (map_it2)
1158 ! update KS matrix
1159 CALL para_env%sum(hmat)
1160 CALL dbcsr_get_block_p(matrix=ksmat, row=atom_a, col=atom_a, block=h_block, found=found)
1161 IF (found) THEN
1162 h_block(1:nsgfa, 1:nsgfa) = h_block(1:nsgfa, 1:nsgfa) + hmat(1:nsgfa, 1:nsgfa)
1163 END IF
1164 IF (calculate_forces) THEN
1165 force(ikind)%rho_elec(:, iatom) = force(ikind)%rho_elec(:, iatom) + 2.0_dp*force_a(:)
1166 IF (use_virial) THEN
1167 IF (use_virial .AND. calculate_forces) THEN
1168 virial%pv_lrigpw = virial%pv_lrigpw + 2.0_dp*my_virial_a
1169 virial%pv_virial = virial%pv_virial + 2.0_dp*my_virial_a
1170 END IF
1171 END IF
1172 END IF
1173 END DO
1174 DEALLOCATE (hab, work, hmat)
1175 IF (calculate_forces) DEALLOCATE (pab, pblk)
1176 END DO
1177
1178 CALL timestop(handle)
1179
1180 END SUBROUTINE integrate_v_rspace_diagonal
1181
1182! **************************************************************************************************
1183!> \brief computes integrals of product of v_rspace times a basis function (vector function)
1184!> and possible forces
1185!> \param qs_env ...
1186!> \param v_rspace ...
1187!> \param f_coef ...
1188!> \param f_integral ...
1189!> \param calculate_forces ...
1190!> \param basis_type ...
1191!> \author JGH [8.2024]
1192! **************************************************************************************************
1193 SUBROUTINE integrate_function(qs_env, v_rspace, f_coef, f_integral, calculate_forces, basis_type)
1194 TYPE(qs_environment_type), POINTER :: qs_env
1195 TYPE(pw_r3d_rs_type), INTENT(IN) :: v_rspace
1196 REAL(kind=dp), DIMENSION(:), INTENT(IN) :: f_coef
1197 REAL(kind=dp), DIMENSION(:), INTENT(OUT) :: f_integral
1198 LOGICAL, INTENT(IN) :: calculate_forces
1199 CHARACTER(len=*), INTENT(IN) :: basis_type
1200
1201 CHARACTER(len=*), PARAMETER :: routinen = 'integrate_function'
1202
1203 INTEGER :: atom_a, group_size, handle, i, iatom, igrid_level, ikind, ipgf, iset, maxco, &
1204 maxsgf_set, my_pos, na1, natom, ncoa, nkind, nseta, offset, sgfa
1205 INTEGER, ALLOCATABLE, DIMENSION(:) :: atom_of_kind
1206 INTEGER, DIMENSION(:), POINTER :: la_max, la_min, npgfa, nsgfa
1207 INTEGER, DIMENSION(:, :), POINTER :: first_sgfa
1208 LOGICAL :: use_virial
1209 REAL(kind=dp) :: eps_rho_rspace, radius
1210 REAL(kind=dp), DIMENSION(3) :: force_a, force_b, ra
1211 REAL(kind=dp), DIMENSION(3, 3) :: my_virial_a, my_virial_b
1212 REAL(kind=dp), DIMENSION(:, :), POINTER :: hab, pab, sphi_a, work, zeta
1213 TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
1214 TYPE(cell_type), POINTER :: cell
1215 TYPE(dft_control_type), POINTER :: dft_control
1216 TYPE(gridlevel_info_type), POINTER :: gridlevel_info
1217 TYPE(gto_basis_set_type), POINTER :: orb_basis_set
1218 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
1219 TYPE(pw_env_type), POINTER :: pw_env
1220 TYPE(qs_force_type), DIMENSION(:), POINTER :: force
1221 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
1222 TYPE(realspace_grid_type), DIMENSION(:), POINTER :: rs_v
1223 TYPE(realspace_grid_type), POINTER :: rs_grid
1224 TYPE(virial_type), POINTER :: virial
1225
1226 CALL timeset(routinen, handle)
1227
1228 CALL get_qs_env(qs_env=qs_env, pw_env=pw_env)
1229 gridlevel_info => pw_env%gridlevel_info
1230
1231 CALL pw_env_get(pw_env, rs_grids=rs_v)
1232 CALL potential_pw2rs(rs_v, v_rspace, pw_env)
1233
1234 CALL get_qs_env(qs_env=qs_env, &
1235 atomic_kind_set=atomic_kind_set, &
1236 qs_kind_set=qs_kind_set, &
1237 force=force, &
1238 cell=cell, &
1239 dft_control=dft_control, &
1240 nkind=nkind, &
1241 natom=natom, &
1242 particle_set=particle_set, &
1243 pw_env=pw_env, &
1244 virial=virial)
1245 CALL get_atomic_kind_set(atomic_kind_set=atomic_kind_set, atom_of_kind=atom_of_kind)
1246
1247 use_virial = virial%pv_availability .AND. (.NOT. virial%pv_numer)
1248 IF (use_virial) THEN
1249 cpabort("Virial NYA")
1250 END IF
1251
1252 eps_rho_rspace = dft_control%qs_control%eps_rho_rspace
1253
1254 CALL get_qs_kind_set(qs_kind_set, &
1255 maxco=maxco, maxsgf_set=maxsgf_set, basis_type=basis_type)
1256 ALLOCATE (hab(maxco, 1), pab(maxco, 1), work(maxco, 1))
1257
1258 offset = 0
1259 my_pos = v_rspace%pw_grid%para%group%mepos
1260 group_size = v_rspace%pw_grid%para%group%num_pe
1261
1262 DO iatom = 1, natom
1263 ikind = particle_set(iatom)%atomic_kind%kind_number
1264 atom_a = atom_of_kind(iatom)
1265 CALL get_qs_kind(qs_kind_set(ikind), basis_set=orb_basis_set, basis_type=basis_type)
1266 CALL get_gto_basis_set(gto_basis_set=orb_basis_set, &
1267 first_sgf=first_sgfa, &
1268 lmax=la_max, &
1269 lmin=la_min, &
1270 npgf=npgfa, &
1271 nset=nseta, &
1272 nsgf_set=nsgfa, &
1273 sphi=sphi_a, &
1274 zet=zeta)
1275 ra(:) = pbc(particle_set(iatom)%r, cell)
1276
1277 force_a(:) = 0._dp
1278 force_b(:) = 0._dp
1279 my_virial_a(:, :) = 0._dp
1280 my_virial_b(:, :) = 0._dp
1281
1282 DO iset = 1, nseta
1283
1284 ncoa = npgfa(iset)*ncoset(la_max(iset))
1285 sgfa = first_sgfa(1, iset)
1286
1287 hab = 0._dp
1288 pab = 0._dp
1289
1290 DO i = 1, nsgfa(iset)
1291 work(i, 1) = f_coef(offset + i)
1292 END DO
1293
1294 CALL dgemm("N", "N", ncoa, 1, nsgfa(iset), &
1295 1.0_dp, sphi_a(1, sgfa), SIZE(sphi_a, 1), &
1296 work(1, 1), SIZE(work, 1), &
1297 0.0_dp, pab(1, 1), SIZE(pab, 1))
1298
1299 DO ipgf = 1, npgfa(iset)
1300
1301 na1 = (ipgf - 1)*ncoset(la_max(iset))
1302
1303 igrid_level = gaussian_gridlevel(gridlevel_info, zeta(ipgf, iset))
1304 rs_grid => rs_v(igrid_level)
1305
1306 IF (map_gaussian_here(rs_grid, cell%h_inv, ra, offset, group_size, my_pos)) THEN
1307 radius = exp_radius_very_extended(la_min=la_min(iset), la_max=la_max(iset), &
1308 lb_min=0, lb_max=0, ra=ra, rb=ra, rp=ra, &
1309 zetp=zeta(ipgf, iset), eps=eps_rho_rspace, &
1310 prefactor=1.0_dp, cutoff=1.0_dp)
1311
1312 IF (calculate_forces) THEN
1313 CALL integrate_pgf_product(la_max=la_max(iset), &
1314 zeta=zeta(ipgf, iset), la_min=la_min(iset), &
1315 lb_max=0, zetb=0.0_dp, lb_min=0, &
1316 ra=ra, rab=(/0.0_dp, 0.0_dp, 0.0_dp/), &
1317 rsgrid=rs_grid, &
1318 hab=hab, pab=pab, o1=na1, o2=0, radius=radius, &
1319 calculate_forces=.true., &
1320 force_a=force_a, force_b=force_b, &
1321 use_virial=use_virial, &
1322 my_virial_a=my_virial_a, my_virial_b=my_virial_b)
1323 ELSE
1324 CALL integrate_pgf_product(la_max=la_max(iset), &
1325 zeta=zeta(ipgf, iset), la_min=la_min(iset), &
1326 lb_max=0, zetb=0.0_dp, lb_min=0, &
1327 ra=ra, rab=(/0.0_dp, 0.0_dp, 0.0_dp/), &
1328 rsgrid=rs_grid, &
1329 hab=hab, o1=na1, o2=0, radius=radius, &
1330 calculate_forces=.false.)
1331 END IF
1332
1333 END IF
1334
1335 END DO
1336 !
1337 CALL dgemm("T", "N", nsgfa(iset), 1, ncoa, 1.0_dp, sphi_a(1, sgfa), &
1338 SIZE(sphi_a, 1), hab(1, 1), SIZE(hab, 1), 0.0_dp, work(1, 1), SIZE(work, 1))
1339 DO i = 1, nsgfa(iset)
1340 f_integral(offset + i) = work(i, 1)
1341 END DO
1342
1343 offset = offset + nsgfa(iset)
1344
1345 END DO
1346
1347 IF (calculate_forces) THEN
1348 force(ikind)%rho_elec(:, atom_a) = force(ikind)%rho_elec(:, atom_a) + force_a(:)
1349 IF (use_virial) THEN
1350 virial%pv_virial = virial%pv_virial + my_virial_a
1351 END IF
1352 END IF
1353
1354 END DO
1355
1356 DEALLOCATE (hab, pab, work)
1357
1358 CALL timestop(handle)
1359
1360 END SUBROUTINE integrate_function
1361
static GRID_HOST_DEVICE int modulo(int a, int m)
Equivalent of Fortran's MODULO, which always return a positive number. https://gcc....
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.
All kind of helpful little routines.
Definition ao_util.F:14
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
Define the atomic kind types and their sub types.
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.
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.
Holds information on atomic properties.
subroutine, public atprop_array_init(atarray, natom)
...
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)
...
Handles all functions related to the CELL.
Definition cell_types.F:15
Defines control structures, which contain the parameters and the settings for the DFT-based calculati...
subroutine, public dbcsr_get_block_p(matrix, row, col, block, found, row_size, col_size)
...
Definition of the atomic potential types.
integer function, public gaussian_gridlevel(gridlevel_info, exponent)
...
Fortran API for the grid package, which is written in C.
Definition grid_api.F:12
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:279
Defines the basic variable types.
Definition kinds.F:23
integer, parameter, public dp
Definition kinds.F:34
contains the types and subroutines for dealing with the lri_env lri : local resolution of the identit...
Utility routines for the memory handling.
Interface to the message passing library MPI.
Provides Cartesian and spherical orbital pointers and indices.
integer, dimension(:), allocatable, public ncoset
integer, dimension(:, :, :), allocatable, public coset
Define the data structure for the particle information.
container for various plainwaves related things
subroutine, public pw_env_get(pw_env, pw_pools, cube_info, gridlevel_info, auxbas_pw_pool, auxbas_grid, auxbas_rs_desc, auxbas_rs_grid, rs_descs, rs_grids, xc_pw_pool, vdw_pw_pool, poisson_env, interp_section)
returns the various attributes of the pw env
subroutine, public get_qs_env(qs_env, atomic_kind_set, qs_kind_set, cell, super_cell, cell_ref, use_ref_cell, kpoints, dft_control, mos, sab_orb, sab_all, qmmm, qmmm_periodic, sac_ae, sac_ppl, sac_lri, sap_ppnl, sab_vdw, sab_scp, sap_oce, sab_lrc, sab_se, sab_xtbe, sab_tbe, sab_core, sab_xb, sab_xtb_pp, sab_xtb_nonbond, sab_almo, sab_kp, sab_kp_nosym, sab_cneo, particle_set, energy, force, matrix_h, matrix_h_im, matrix_ks, matrix_ks_im, matrix_vxc, run_rtp, rtp, matrix_h_kp, matrix_h_im_kp, matrix_ks_kp, matrix_ks_im_kp, matrix_vxc_kp, kinetic_kp, matrix_s_kp, matrix_w_kp, matrix_s_ri_aux_kp, matrix_s, matrix_s_ri_aux, matrix_w, matrix_p_mp2, matrix_p_mp2_admm, rho, rho_xc, pw_env, ewald_env, ewald_pw, active_space, mpools, input, para_env, blacs_env, scf_control, rel_control, kinetic, qs_charges, vppl, rho_core, rho_nlcc, rho_nlcc_g, ks_env, ks_qmmm_env, wf_history, scf_env, local_particles, local_molecules, distribution_2d, dbcsr_dist, molecule_kind_set, molecule_set, subsys, cp_subsys, oce, local_rho_set, rho_atom_set, task_list, task_list_soft, rho0_atom_set, rho0_mpole, rhoz_set, rhoz_cneo_set, ecoul_1c, rho0_s_rs, rho0_s_gs, rhoz_cneo_s_rs, rhoz_cneo_s_gs, do_kpoints, has_unit_metric, requires_mo_derivs, mo_derivs, mo_loc_history, nkind, natom, nelectron_total, nelectron_spin, efield, neighbor_list_id, linres_control, xas_env, virial, cp_ddapc_env, cp_ddapc_ewald, outer_scf_history, outer_scf_ihistory, x_data, et_coupling, dftb_potential, results, se_taper, se_store_int_env, se_nddo_mpole, se_nonbond_env, admm_env, lri_env, lri_density, exstate_env, ec_env, harris_env, dispersion_env, gcp_env, vee, rho_external, external_vxc, mask, mp2_env, bs_env, kg_env, wanniercentres, atprop, ls_scf_env, do_transport, transport_env, v_hartree_rspace, s_mstruct_changed, rho_changed, potential_changed, forces_up_to_date, mscfg_env, almo_scf_env, gradient_history, variable_history, embed_pot, spin_embed_pot, polar_env, mos_last_converged, eeq, rhs, do_rixs, tb_tblite)
Get the QUICKSTEP environment.
Build up the plane wave density by collocating the primitive Gaussian functions (pgf).
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...
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
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...
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
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
subroutine, public integrate_ppl_rspace(rho_rspace, qs_env)
computes the forces/virial due to the local pseudopotential
subroutine, public integrate_rho_nlcc(rho_rspace, qs_env)
computes the forces/virial due to the nlcc pseudopotential
Define the quickstep kind type and their sub types.
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, floating, name, element_symbol, pao_basis_size, pao_model_file, pao_potentials, pao_descriptors, nelec)
Get attributes of an atomic kind.
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.
subroutine, public transfer_pw2rs(rs, pw)
...
pure logical function, public map_gaussian_here(rs_grid, h_inv, ra, offset, group_size, my_pos)
...
subroutine, public rs_grid_zero(rs)
Initialize grid to zero.
Transfers densities from PW to RS grids and potentials from PW to RS.
subroutine, public potential_pw2rs(rs_v, v_rspace, pw_env)
transfers a potential from a pw_grid to a vector of realspace multigrids
Provides all information about an atomic kind.
type for the atomic properties
Type defining parameters related to the simulation cell.
Definition cell_types.F:55
stores all the informations relevant to an mpi environment
contained for different pw related things
Provides all information about a quickstep kind.