(git:d1312bc)
Loading...
Searching...
No Matches
hartree_local_methods.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! **************************************************************************************************
13 USE cell_types, ONLY: cell_type
19 USE kinds, ONLY: dp
20 USE mathconstants, ONLY: fourpi,&
21 pi
23 USE orbital_pointers, ONLY: indso,&
24 nsoset
25 USE pw_env_types, ONLY: pw_env_get,&
39 USE qs_kind_types, ONLY: get_qs_kind,&
45 USE qs_rho0_types, ONLY: get_rho0_mpole,&
51 USE util, ONLY: get_limit
52#include "./base/base_uses.f90"
53
54 IMPLICIT NONE
55
56 PRIVATE
57
58 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'hartree_local_methods'
59
60 ! Public Subroutine
61
63
64CONTAINS
65
66! **************************************************************************************************
67!> \brief ...
68!> \param hartree_local ...
69!> \param natom ...
70! **************************************************************************************************
71 SUBROUTINE init_coulomb_local(hartree_local, natom)
72
73 TYPE(hartree_local_type), POINTER :: hartree_local
74 INTEGER, INTENT(IN) :: natom
75
76 CHARACTER(len=*), PARAMETER :: routinen = 'init_coulomb_local'
77
78 INTEGER :: handle
79 TYPE(ecoul_1center_type), DIMENSION(:), POINTER :: ecoul_1c
80
81 CALL timeset(routinen, handle)
82
83 NULLIFY (ecoul_1c)
84 ! Allocate and Initialize 1-center Potentials and Integrals
85 CALL allocate_ecoul_1center(ecoul_1c, natom)
86 hartree_local%ecoul_1c => ecoul_1c
87
88 CALL timestop(handle)
89
90 END SUBROUTINE init_coulomb_local
91
92! **************************************************************************************************
93!> \brief Calculates Hartree potential for hard and soft densities (including
94!> nuclear charge and compensation charges) using numerical integration
95!> \param vrad_h ...
96!> \param vrad_s ...
97!> \param rrad_h ...
98!> \param rrad_s ...
99!> \param rrad_0 ...
100!> \param rrad_z ...
101!> \param grid_atom ...
102!> \par History
103!> 05.2012 JGH refactoring
104!> \author ??
105! **************************************************************************************************
106 SUBROUTINE calculate_vh_1center(vrad_h, vrad_s, rrad_h, rrad_s, rrad_0, rrad_z, grid_atom)
107
108 REAL(dp), DIMENSION(:, :), INTENT(INOUT) :: vrad_h, vrad_s
109 TYPE(rho_atom_coeff), DIMENSION(:), INTENT(IN) :: rrad_h, rrad_s
110 REAL(dp), DIMENSION(:, :), INTENT(IN) :: rrad_0
111 REAL(dp), DIMENSION(:), INTENT(IN) :: rrad_z
112 TYPE(grid_atom_type), POINTER :: grid_atom
113
114 CHARACTER(len=*), PARAMETER :: routinen = 'calculate_Vh_1center'
115
116 INTEGER :: handle, ir, iso, ispin, l_ang, &
117 max_s_harm, nchannels, nr, nspins
118 REAL(dp) :: i1_down, i1_up, i2_down, i2_up, prefactor
119 REAL(dp), ALLOCATABLE, DIMENSION(:) :: rho_1, rho_2
120 REAL(dp), DIMENSION(:), POINTER :: wr
121 REAL(dp), DIMENSION(:, :), POINTER :: oor2l, r2l
122
123 CALL timeset(routinen, handle)
124
125 nr = grid_atom%nr
126 max_s_harm = SIZE(vrad_h, 2)
127 nspins = SIZE(rrad_h, 1)
128 nchannels = SIZE(rrad_0, 2)
129
130 r2l => grid_atom%rad2l
131 oor2l => grid_atom%oorad2l
132 wr => grid_atom%wr
133
134 ALLOCATE (rho_1(nr), rho_2(nr))
135
136 DO iso = 1, max_s_harm
137 rho_1(:) = 0.0_dp
138 rho_2(:) = 0.0_dp
139 IF (iso == 1) rho_1(:) = rrad_z(:)
140 IF (iso <= nchannels) rho_2(:) = rrad_0(:, iso)
141 DO ispin = 1, nspins
142 rho_1(:) = rho_1(:) + rrad_h(ispin)%r_coef(:, iso)
143 rho_2(:) = rho_2(:) + rrad_s(ispin)%r_coef(:, iso)
144 END DO
145
146 l_ang = indso(1, iso)
147 prefactor = fourpi/(2._dp*l_ang + 1._dp)
148
149 rho_1(:) = rho_1(:)*wr(:)
150 rho_2(:) = rho_2(:)*wr(:)
151
152 i1_up = 0.0_dp
153 i1_down = 0.0_dp
154 i2_up = 0.0_dp
155 i2_down = 0.0_dp
156
157 i1_up = r2l(nr, l_ang)*rho_1(nr)
158 i2_up = r2l(nr, l_ang)*rho_2(nr)
159
160 DO ir = nr - 1, 1, -1
161 i1_down = i1_down + oor2l(ir, l_ang + 1)*rho_1(ir)
162 i2_down = i2_down + oor2l(ir, l_ang + 1)*rho_2(ir)
163 END DO
164
165 vrad_h(nr, iso) = vrad_h(nr, iso) + prefactor* &
166 (oor2l(nr, l_ang + 1)*i1_up + r2l(nr, l_ang)*i1_down)
167 vrad_s(nr, iso) = vrad_s(nr, iso) + prefactor* &
168 (oor2l(nr, l_ang + 1)*i2_up + r2l(nr, l_ang)*i2_down)
169
170 DO ir = nr - 1, 1, -1
171 i1_up = i1_up + r2l(ir, l_ang)*rho_1(ir)
172 i1_down = i1_down - oor2l(ir, l_ang + 1)*rho_1(ir)
173 i2_up = i2_up + r2l(ir, l_ang)*rho_2(ir)
174 i2_down = i2_down - oor2l(ir, l_ang + 1)*rho_2(ir)
175
176 vrad_h(ir, iso) = vrad_h(ir, iso) + prefactor* &
177 (oor2l(ir, l_ang + 1)*i1_up + r2l(ir, l_ang)*i1_down)
178 vrad_s(ir, iso) = vrad_s(ir, iso) + prefactor* &
179 (oor2l(ir, l_ang + 1)*i2_up + r2l(ir, l_ang)*i2_down)
180
181 END DO
182
183 END DO
184
185 DEALLOCATE (rho_1, rho_2)
186
187 CALL timestop(handle)
188
189 END SUBROUTINE calculate_vh_1center
190
191! **************************************************************************************************
192!> \brief Calculates one center GAPW Hartree energies and matrix elements
193!> Hartree potentials are input
194!> Takes possible background charge into account
195!> Special case for densities without core charge
196!> \param qs_env ...
197!> \param energy_hartree_1c ...
198!> \param ecoul_1c ...
199!> \param local_rho_set ...
200!> \param para_env ...
201!> \param tddft ...
202!> \param local_rho_set_2nd ...
203!> \param core_2nd ...
204!> \par History
205!> 05.2012 JGH refactoring
206!> \author ??
207! **************************************************************************************************
208 SUBROUTINE vh_1c_gg_integrals(qs_env, energy_hartree_1c, ecoul_1c, local_rho_set, para_env, tddft, local_rho_set_2nd, &
209 core_2nd)
210
211 TYPE(qs_environment_type), POINTER :: qs_env
212 REAL(kind=dp), INTENT(out) :: energy_hartree_1c
213 TYPE(ecoul_1center_type), DIMENSION(:), POINTER :: ecoul_1c
214 TYPE(local_rho_type), POINTER :: local_rho_set
215 TYPE(mp_para_env_type), POINTER :: para_env
216 LOGICAL, INTENT(IN) :: tddft
217 TYPE(local_rho_type), OPTIONAL, POINTER :: local_rho_set_2nd
218 LOGICAL, INTENT(IN), OPTIONAL :: core_2nd
219
220 CHARACTER(LEN=*), PARAMETER :: routinen = 'Vh_1c_gg_integrals'
221
222 INTEGER :: bo(2), handle, iat, iatom, ikind, ipgf1, is1, iset1, iso, l_ang, llmax, &
223 llmax_nuc, lmax0, lmax0_2nd, lmax_0, m1, max_iso, max_iso_not0, max_iso_not0_nuc, &
224 max_s_harm, max_s_harm_nuc, maxl, maxl_nuc, maxso, maxso_nuc, mepos, n1, nat, nchan_0, &
225 nkind, nr, nset, nset_nuc, nsotot, nsotot_nuc, nspins, num_pe
226 INTEGER, ALLOCATABLE, DIMENSION(:) :: cg_n_list
227 INTEGER, ALLOCATABLE, DIMENSION(:, :, :) :: cg_list
228 INTEGER, DIMENSION(:), POINTER :: atom_list, lmax, lmax_nuc, lmin, &
229 lmin_nuc, npgf, npgf_nuc
230 LOGICAL :: cneo, core_charge, l_2nd_local_rho, &
231 my_core_2nd, my_periodic, paw_atom
232 REAL(dp) :: back_ch, ecoul_1_z_cneo, factor
233 REAL(dp), ALLOCATABLE, DIMENSION(:) :: gexp, sqrtwr
234 REAL(dp), ALLOCATABLE, DIMENSION(:, :) :: avh1b_00, avh1b_00_nuc, avh1b_hh, &
235 avh1b_hh_nuc, avh1b_ss, avh1b_ss_nuc, &
236 g0_h_w
237 REAL(dp), DIMENSION(:), POINTER :: rrad_z, vrrad_z
238 REAL(dp), DIMENSION(:, :), POINTER :: g0_h, g0_h_2nd, gsph, gsph_nuc, rrad_0, &
239 vh1_h, vh1_s, vrrad_0, zet, zet_nuc
240 REAL(dp), DIMENSION(:, :, :), POINTER :: my_cg, qlm_gg, qlm_gg_2nd
241 TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
242 TYPE(cell_type), POINTER :: cell
243 TYPE(cneo_potential_type), POINTER :: cneo_potential
244 TYPE(dft_control_type), POINTER :: dft_control
245 TYPE(grid_atom_type), POINTER :: grid_atom
246 TYPE(gto_basis_set_type), POINTER :: basis_1c, nuc_basis
247 TYPE(harmonics_atom_type), POINTER :: harmonics
248 TYPE(pw_env_type), POINTER :: pw_env
249 TYPE(pw_poisson_type), POINTER :: poisson_env
250 TYPE(qs_charges_type), POINTER :: qs_charges
251 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
252 TYPE(rho0_atom_type), DIMENSION(:), POINTER :: rho0_atom_set, rho0_atom_set_2nd
253 TYPE(rho0_mpole_type), POINTER :: rho0_mpole, rho0_mpole_2nd
254 TYPE(rho_atom_type), DIMENSION(:), POINTER :: rho_atom_set, rho_atom_set_2nd
255 TYPE(rho_atom_type), POINTER :: rho_atom
256 TYPE(rhoz_cneo_type), DIMENSION(:), POINTER :: rhoz_cneo_set
257 TYPE(rhoz_cneo_type), POINTER :: rhoz_cneo
258 TYPE(rhoz_type), DIMENSION(:), POINTER :: rhoz_set, rhoz_set_2nd
259
260 CALL timeset(routinen, handle)
261
262 NULLIFY (cell, dft_control, poisson_env, pw_env, qs_charges)
263 NULLIFY (atomic_kind_set, qs_kind_set, rho_atom_set, rho0_atom_set)
264 NULLIFY (rho0_mpole, rhoz_set)
265 NULLIFY (atom_list, grid_atom, harmonics)
266 NULLIFY (basis_1c, lmin, lmax, npgf, zet)
267 NULLIFY (gsph)
268 NULLIFY (rhoz_cneo, rhoz_cneo_set)
269
270 CALL get_qs_env(qs_env=qs_env, &
271 cell=cell, dft_control=dft_control, &
272 atomic_kind_set=atomic_kind_set, &
273 qs_kind_set=qs_kind_set, &
274 pw_env=pw_env, qs_charges=qs_charges)
275
276 CALL pw_env_get(pw_env, poisson_env=poisson_env)
277 my_periodic = (poisson_env%method == pw_poisson_periodic)
278
279 back_ch = qs_charges%background*cell%deth
280
281 ! rhoz_set is not accessed in TDDFT
282 CALL get_local_rho(local_rho_set, rho_atom_set, rho0_atom_set, rho0_mpole, rhoz_set, &
283 rhoz_cneo_set) ! for integral space
284
285 ! for forces we need a second local_rho_set
286 l_2nd_local_rho = .false.
287 IF (PRESENT(local_rho_set_2nd)) THEN ! for potential
288 l_2nd_local_rho = .true.
289 NULLIFY (rho_atom_set_2nd, rho0_atom_set_2nd, rhoz_set_2nd) ! for potential
290 CALL get_local_rho(local_rho_set_2nd, rho_atom_set_2nd, rho0_atom_set_2nd, rho0_mpole_2nd, rhoz_set=rhoz_set_2nd)
291 END IF
292
293 nkind = SIZE(atomic_kind_set, 1)
294 nspins = dft_control%nspins
295
296 core_charge = .NOT. tddft ! for forces mixed version
297 my_core_2nd = .true.
298 IF (PRESENT(core_2nd)) my_core_2nd = .NOT. core_2nd ! if my_core_2nd true, include core charge
299
300 ! The aim of the following code was to return immediately if the subroutine
301 ! was called for triplet excited states in spin-restricted case. This check
302 ! is also performed before invocation of this subroutine. It should be save
303 ! to remove the optional argument 'do_triplet' from the subroutine interface.
304 !IF (tddft) THEN
305 ! CPASSERT(PRESENT(do_triplet))
306 ! IF (nspins == 1 .AND. do_triplet) RETURN
307 !END IF
308
309 CALL get_qs_kind_set(qs_kind_set, maxg_iso_not0=max_iso)
310 CALL get_rho0_mpole(rho0_mpole=rho0_mpole, lmax_0=lmax_0)
311
312 ! Put to 0 the local hartree energy contribution from 1 center integrals
313 energy_hartree_1c = 0.0_dp
314 ! Restore total quantum nuclear density to zero
315 qs_charges%total_rho1_hard_nuc = 0.0_dp
316 rho0_mpole%tot_rhoz_cneo_s = 0.0_dp
317
318 ! Here starts the loop over all the atoms
319 DO ikind = 1, nkind
320
321 CALL get_atomic_kind(atomic_kind_set(ikind), atom_list=atom_list, natom=nat)
322 NULLIFY (cneo_potential)
323 CALL get_qs_kind(qs_kind_set(ikind), &
324 grid_atom=grid_atom, &
325 harmonics=harmonics, ngrid_rad=nr, &
326 max_iso_not0=max_iso_not0, paw_atom=paw_atom, &
327 cneo_potential=cneo_potential)
328 CALL get_qs_kind(qs_kind_set(ikind), &
329 basis_set=basis_1c, basis_type="GAPW_1C")
330
331 cneo = ASSOCIATED(cneo_potential)
332 IF (cneo .AND. tddft) THEN
333 cpabort("Electronic TDDFT with CNEO quantum nuclei is not implemented.")
334 END IF
335
336 NULLIFY (nuc_basis)
337 max_iso_not0_nuc = 0
338 IF (cneo) THEN
339 cpassert(paw_atom)
340 CALL get_qs_kind(qs_kind_set(ikind), &
341 basis_set=nuc_basis, basis_type="NUC")
342 max_iso_not0_nuc = cneo_potential%harmonics%max_iso_not0
343 END IF
344
345 IF (paw_atom) THEN
346!=========== PAW ===============
347 CALL get_gto_basis_set(gto_basis_set=basis_1c, lmax=lmax, lmin=lmin, &
348 maxso=maxso, npgf=npgf, maxl=maxl, &
349 nset=nset, zet=zet)
350
351 max_s_harm = harmonics%max_s_harm
352 llmax = harmonics%llmax
353
354 nsotot = maxso*nset
355 ALLOCATE (gsph(nr, nsotot))
356 ALLOCATE (gexp(nr))
357 ALLOCATE (sqrtwr(nr), g0_h_w(nr, 0:lmax_0))
358
359 ALLOCATE (avh1b_hh(nsotot, nsotot))
360 ALLOCATE (avh1b_ss(nsotot, nsotot))
361 ALLOCATE (avh1b_00(nsotot, nsotot))
362
363 NULLIFY (qlm_gg, g0_h)
364 CALL get_rho0_mpole(rho0_mpole=rho0_mpole, ikind=ikind, &
365 l0_ikind=lmax0, &
366 qlm_gg=qlm_gg, g0_h=g0_h) ! Qlm_gg of density
367
368 IF (PRESENT(local_rho_set_2nd)) THEN ! for potential
369 NULLIFY (qlm_gg_2nd, g0_h_2nd)
370 CALL get_rho0_mpole(rho0_mpole=rho0_mpole_2nd, ikind=ikind, &
371 l0_ikind=lmax0_2nd, &
372 qlm_gg=qlm_gg_2nd, g0_h=g0_h_2nd) ! Qlm_gg of density
373 END IF
374 nchan_0 = nsoset(lmax0)
375
376 IF (nchan_0 > max(max_iso_not0, max_iso_not0_nuc)) THEN
377 cpabort("channels for rho0 > # max of spherical harmonics")
378 END IF
379
380 NULLIFY (vh1_h, vh1_s)
381 ALLOCATE (vh1_h(nr, max_iso_not0))
382 ALLOCATE (vh1_s(nr, max(max_iso_not0, max_iso_not0_nuc, nchan_0)))
383
384 NULLIFY (lmax_nuc, lmin_nuc, npgf_nuc, zet_nuc, gsph_nuc)
385 maxso_nuc = 0
386 maxl_nuc = -1
387 nset_nuc = 0
388 max_s_harm_nuc = 0
389 llmax_nuc = -1
390 nsotot_nuc = 0
391 IF (cneo) THEN
392 CALL get_gto_basis_set(gto_basis_set=nuc_basis, lmax=lmax_nuc, &
393 lmin=lmin_nuc, maxso=maxso_nuc, npgf=npgf_nuc, &
394 maxl=maxl_nuc, nset=nset_nuc, zet=zet_nuc)
395
396 max_s_harm_nuc = cneo_potential%harmonics%max_s_harm
397 llmax_nuc = cneo_potential%harmonics%llmax
398 nsotot_nuc = maxso_nuc*nset_nuc
399 ALLOCATE (gsph_nuc(nr, nsotot_nuc))
400 gsph_nuc = 0.0_dp
401 ALLOCATE (avh1b_hh_nuc(nsotot_nuc, nsotot_nuc))
402 ALLOCATE (avh1b_ss_nuc(nsotot_nuc, nsotot_nuc))
403 ALLOCATE (avh1b_00_nuc(nsotot_nuc, nsotot_nuc))
404 END IF
405
406 ALLOCATE (cg_list(2, nsoset(max(maxl, maxl_nuc))**2, &
407 max(max_s_harm, max_s_harm_nuc)), &
408 cg_n_list(max(max_s_harm, max_s_harm_nuc)))
409
410 NULLIFY (rrad_z, my_cg)
411 my_cg => harmonics%my_CG
412
413 ! set to zero temporary arrays
414 sqrtwr = 0.0_dp
415 g0_h_w = 0.0_dp
416 gexp = 0.0_dp
417 gsph = 0.0_dp
418
419 sqrtwr(1:nr) = sqrt(grid_atom%wr(1:nr))
420 DO l_ang = 0, lmax0
421 g0_h_w(1:nr, l_ang) = g0_h(1:nr, l_ang)*grid_atom%wr(1:nr)
422 END DO
423
424 m1 = 0
425 DO iset1 = 1, nset
426 n1 = nsoset(lmax(iset1))
427 DO ipgf1 = 1, npgf(iset1)
428 gexp(1:nr) = exp(-zet(ipgf1, iset1)*grid_atom%rad2(1:nr))*sqrtwr(1:nr)
429 DO is1 = nsoset(lmin(iset1) - 1) + 1, nsoset(lmax(iset1))
430 iso = is1 + (ipgf1 - 1)*n1 + m1
431 l_ang = indso(1, is1)
432 gsph(1:nr, iso) = grid_atom%rad2l(1:nr, l_ang)*gexp(1:nr)
433 END DO ! is1
434 END DO ! ipgf1
435 m1 = m1 + maxso
436 END DO ! iset1
437
438 IF (cneo) THEN
439 ! initialize nuclear pmat, cpc and e_core to zero
440 DO iat = 1, nat
441 iatom = atom_list(iat)
442 rhoz_cneo_set(iatom)%pmat = 0.0_dp
443 rhoz_cneo_set(iatom)%cpc_h = 0.0_dp
444 rhoz_cneo_set(iatom)%cpc_s = 0.0_dp
445 rhoz_cneo_set(iatom)%e_core = 0.0_dp
446 rhoz_cneo_set(iatom)%ready = .false.
447 END DO
448
449 ! calculate nuclear gsph
450 m1 = 0
451 DO iset1 = 1, nset_nuc
452 n1 = nsoset(lmax_nuc(iset1))
453 DO ipgf1 = 1, npgf_nuc(iset1)
454 gexp(1:nr) = exp(-zet_nuc(ipgf1, iset1)*grid_atom%rad2(1:nr))*sqrtwr(1:nr)
455 DO is1 = nsoset(lmin_nuc(iset1) - 1) + 1, nsoset(lmax_nuc(iset1))
456 iso = is1 + (ipgf1 - 1)*n1 + m1
457 l_ang = indso(1, is1)
458 gsph_nuc(1:nr, iso) = cneo_potential%rad2l(1:nr, l_ang)*gexp(1:nr)
459 END DO ! is1
460 END DO ! ipgf1
461 m1 = m1 + maxso_nuc
462 END DO ! iset1
463 END IF
464
465 ! Distribute the atoms of this kind
466 num_pe = para_env%num_pe
467 mepos = para_env%mepos
468 bo = get_limit(nat, num_pe, mepos)
469
470 DO iat = bo(1), bo(2) !1,nat
471 iatom = atom_list(iat)
472 rho_atom => rho_atom_set(iatom)
473
474 NULLIFY (rrad_z, vrrad_z, rrad_0, vrrad_0)
475 IF (core_charge .AND. .NOT. cneo) THEN
476 rrad_z => rhoz_set(ikind)%r_coef ! for density
477 END IF
478 IF (my_core_2nd .AND. .NOT. cneo) THEN
479 IF (l_2nd_local_rho) THEN
480 vrrad_z => rhoz_set_2nd(ikind)%vr_coef ! for potential
481 ELSE
482 vrrad_z => rhoz_set(ikind)%vr_coef ! for potential
483 END IF
484 END IF
485 rrad_0 => rho0_atom_set(iatom)%rho0_rad_h%r_coef ! for density
486 vrrad_0 => rho0_atom_set(iatom)%vrho0_rad_h%r_coef
487 IF (l_2nd_local_rho) THEN
488 rho_atom => rho_atom_set_2nd(iatom)
489 vrrad_0 => rho0_atom_set_2nd(iatom)%vrho0_rad_h%r_coef ! for potential
490 END IF
491 IF (my_periodic .AND. back_ch > 1.e-3_dp) THEN
492 factor = -2.0_dp*pi/3.0_dp*sqrt(fourpi)*qs_charges%background
493 ELSE
494 factor = 0._dp
495 END IF
496
497 CALL vh_1c_atom_potential(rho_atom, vrrad_0, &
498 grid_atom, my_core_2nd .AND. .NOT. cneo, & ! core charge for potential (2nd)
499 vrrad_z, vh1_h, vh1_s, &
500 nchan_0, nspins, max_iso_not0, factor)
501
502 IF (l_2nd_local_rho) rho_atom => rho_atom_set(iatom) ! rho_atom for density
503
504 ecoul_1_z_cneo = 0.0_dp
505 IF (cneo) THEN
506 rhoz_cneo => rhoz_cneo_set(iatom)
507 ! Add the soft tail of nuclear Hartree potential to total Vh1_s first.
508 ! vrho already contains the -zeff factor.
509 DO iso = 1, max_iso_not0_nuc
510 vh1_s(:, iso) = vh1_s(:, iso) + rhoz_cneo%vrho_rad_s(:, iso)
511 END DO
512
513 ! Build nuclear 1c integrals according to Vh1_h from electronic density only
514 ! and Vh1_s from electron density and last step (or initial guess) nuclear density
515 ! Vh1_h_nuc = -Z*Vh1_h, Vh1_s_nuc = -Z*Vh1_s
516 CALL vh_1c_nuc_integrals(rhoz_cneo, cneo_potential%zeff, &
517 avh1b_hh_nuc, avh1b_ss_nuc, avh1b_00_nuc, vh1_h, vh1_s, &
518 max_iso_not0, max_iso_not0_nuc, &
519 max_s_harm_nuc, llmax_nuc, cg_list, cg_n_list, &
520 nset_nuc, npgf_nuc, lmin_nuc, lmax_nuc, nsotot_nuc, maxso_nuc, &
521 nchan_0, gsph_nuc, g0_h_w, cneo_potential%harmonics%my_CG, &
522 cneo_potential%Qlm_gg)
523
524 ! Solve the nuclear 1c problem
525 CALL calculate_rhoz_cneo(rhoz_cneo, cneo_potential, cg_list, cg_n_list, nset_nuc, &
526 npgf_nuc, lmin_nuc, lmax_nuc, maxl_nuc, maxso_nuc)
527 ! slm_int(iso=1) = sqrt(4*Pi), slm_int(iso>1) = 0, without using Lebedev grid
528 ! when printing, nuclear density is positive, thus needing the minus sign
529 qs_charges%total_rho1_hard_nuc = qs_charges%total_rho1_hard_nuc - sqrt(fourpi) &
530 *sum(rhoz_cneo%rho_rad_h(:, 1)*grid_atom%wr(:))
531 rho0_mpole%tot_rhoz_cneo_s = rho0_mpole%tot_rhoz_cneo_s - sqrt(fourpi) &
532 *sum(rhoz_cneo%rho_rad_s(:, 1)*grid_atom%wr(:))
533
534 ! Calculate the contributions to Ecoul coming from Vh1_h*rhoz_cneo_h and Vh1_s*rhoz_cneo_s.
535 ! Self-interaction of the quantum nucleus is already removed in ecoul_1_z_cneo.
536 ! rho already contains the -zeff factor.
537 DO iso = 1, min(max_iso_not0, max_iso_not0_nuc)
538 ecoul_1_z_cneo = ecoul_1_z_cneo + 0.5_dp* &
539 (sum(vh1_h(:, iso)*rhoz_cneo%rho_rad_h(:, iso)*grid_atom%wr(:)) &
540 - sum(vh1_s(:, iso)*rhoz_cneo%rho_rad_s(:, iso)*grid_atom%wr(:)))
541 END DO
542 DO iso = max_iso_not0 + 1, max_iso_not0_nuc
543 ecoul_1_z_cneo = ecoul_1_z_cneo - 0.5_dp* &
544 sum(vh1_s(:, iso)*rhoz_cneo%rho_rad_s(:, iso)*grid_atom%wr(:))
545 END DO
546
547 ! Add nuclear Hartree potential to total Vh1_h after solving the nuclear 1c problem
548 ! to avoid nuclear self-interaction.
549 ! vrho already contains the -zeff factor.
550 ! Here the min of two max_iso_not0's is chosen, because even when the nuclear one
551 ! is larger, it is meaningless to let Vh1_h have higher angular momentum components
552 ! as Vh1_h now is only used for the electronic part.
553 DO iso = 1, min(max_iso_not0, max_iso_not0_nuc)
554 vh1_h(:, iso) = vh1_h(:, iso) + rhoz_cneo%vrho_rad_h(:, iso)
555 END DO
556 END IF
557
558 CALL vh_1c_atom_energy(energy_hartree_1c, ecoul_1c, rho_atom, rrad_0, &
559 grid_atom, iatom, core_charge .AND. .NOT. cneo, & ! core charge for density
560 rrad_z, vh1_h, vh1_s, nchan_0, nspins, max_iso_not0)
561
562 IF (l_2nd_local_rho) rho_atom => rho_atom_set_2nd(iatom) ! rho_atom for potential (2nd)
563
564 IF (cneo) THEN
565 CALL set_ecoul_1c(ecoul_1c, iatom, ecoul_1_z=ecoul_1c(iatom)%ecoul_1_z + ecoul_1_z_cneo)
566 energy_hartree_1c = energy_hartree_1c + ecoul_1_z_cneo
567 END IF
568
569 CALL vh_1c_atom_integrals(rho_atom, & ! results (int_local_h and int_local_s) written on rho_atom_2nd
570 ! int_local_h and int_local_s are used in update_ks_atom
571 ! on int_local_h mixed core / non-core
572 avh1b_hh, avh1b_ss, avh1b_00, vh1_h, vh1_s, max_iso_not0, &
573 max_s_harm, llmax, cg_list, cg_n_list, &
574 nset, npgf, lmin, lmax, nsotot, maxso, nspins, nchan_0, gsph, &
575 g0_h_w, my_cg, qlm_gg) ! Qlm_gg for density from local_rho_set
576
577 END DO ! iat
578
579 DEALLOCATE (avh1b_hh)
580 DEALLOCATE (avh1b_ss)
581 DEALLOCATE (avh1b_00)
582 IF (ALLOCATED(avh1b_hh_nuc)) DEALLOCATE (avh1b_hh_nuc)
583 IF (ALLOCATED(avh1b_ss_nuc)) DEALLOCATE (avh1b_ss_nuc)
584 IF (ALLOCATED(avh1b_00_nuc)) DEALLOCATE (avh1b_00_nuc)
585 DEALLOCATE (vh1_h, vh1_s)
586 DEALLOCATE (cg_list, cg_n_list)
587 DEALLOCATE (gsph)
588 IF (ASSOCIATED(gsph_nuc)) DEALLOCATE (gsph_nuc)
589 DEALLOCATE (gexp)
590 DEALLOCATE (sqrtwr, g0_h_w)
591
592 IF (cneo) THEN
593 ! broadcast nuclear pmat, cpc and e_core
594 DO iat = 1, nat
595 iatom = atom_list(iat)
596 CALL para_env%sum(rhoz_cneo_set(iatom)%pmat)
597 CALL para_env%sum(rhoz_cneo_set(iatom)%cpc_h)
598 CALL para_env%sum(rhoz_cneo_set(iatom)%cpc_s)
599 CALL para_env%sum(rhoz_cneo_set(iatom)%e_core)
600 rhoz_cneo_set(iatom)%ready = .true.
601 END DO
602 END IF
603 ELSE
604!=========== NO PAW ===============
605! This term is taken care of using the core density as in GPW
606 cycle
607 END IF ! paw
608 END DO ! ikind
609
610 CALL para_env%sum(energy_hartree_1c)
611 CALL para_env%sum(qs_charges%total_rho1_hard_nuc)
612 CALL para_env%sum(rho0_mpole%tot_rhoz_cneo_s)
613
614 CALL timestop(handle)
615
616 END SUBROUTINE vh_1c_gg_integrals
617
618! **************************************************************************************************
619
620! **************************************************************************************************
621!> \brief ...
622!> \param rho_atom ...
623!> \param vrrad_0 ...
624!> \param grid_atom ...
625!> \param core_charge ...
626!> \param vrrad_z ...
627!> \param Vh1_h ...
628!> \param Vh1_s ...
629!> \param nchan_0 ...
630!> \param nspins ...
631!> \param max_iso_not0 ...
632!> \param bfactor ...
633! **************************************************************************************************
634 SUBROUTINE vh_1c_atom_potential(rho_atom, vrrad_0, &
635 grid_atom, core_charge, vrrad_z, Vh1_h, Vh1_s, &
636 nchan_0, nspins, max_iso_not0, bfactor)
637
638 TYPE(rho_atom_type), POINTER :: rho_atom
639 REAL(dp), DIMENSION(:, :), POINTER :: vrrad_0
640 TYPE(grid_atom_type), POINTER :: grid_atom
641 LOGICAL, INTENT(IN) :: core_charge
642 REAL(dp), DIMENSION(:), POINTER :: vrrad_z
643 REAL(dp), DIMENSION(:, :), POINTER :: vh1_h, vh1_s
644 INTEGER, INTENT(IN) :: nchan_0, nspins, max_iso_not0
645 REAL(dp), INTENT(IN) :: bfactor
646
647 INTEGER :: ir, iso, ispin, nr
648 TYPE(rho_atom_coeff), DIMENSION(:), POINTER :: vr_h, vr_s
649
650 nr = grid_atom%nr
651
652 NULLIFY (vr_h, vr_s)
653 CALL get_rho_atom(rho_atom=rho_atom, vrho_rad_h=vr_h, vrho_rad_s=vr_s)
654
655 vh1_h = 0.0_dp
656 vh1_s = 0.0_dp
657
658 IF (core_charge) vh1_h(:, 1) = vrrad_z(:)
659
660 DO iso = 1, nchan_0
661 vh1_s(:, iso) = vrrad_0(:, iso)
662 END DO
663
664 DO ispin = 1, nspins
665 DO iso = 1, max_iso_not0
666 vh1_h(:, iso) = vh1_h(:, iso) + vr_h(ispin)%r_coef(:, iso)
667 vh1_s(:, iso) = vh1_s(:, iso) + vr_s(ispin)%r_coef(:, iso)
668 END DO
669 END DO
670
671 IF (bfactor /= 0._dp) THEN
672 DO ir = 1, nr
673 vh1_h(ir, 1) = vh1_h(ir, 1) + bfactor*grid_atom%rad2(ir)*grid_atom%wr(ir)
674 vh1_s(ir, 1) = vh1_s(ir, 1) + bfactor*grid_atom%rad2(ir)*grid_atom%wr(ir)
675 END DO
676 END IF
677
678 END SUBROUTINE vh_1c_atom_potential
679
680! **************************************************************************************************
681
682! **************************************************************************************************
683!> \brief ...
684!> \param energy_hartree_1c ...
685!> \param ecoul_1c ...
686!> \param rho_atom ...
687!> \param rrad_0 ...
688!> \param grid_atom ...
689!> \param iatom ...
690!> \param core_charge ...
691!> \param rrad_z ...
692!> \param Vh1_h ...
693!> \param Vh1_s ...
694!> \param nchan_0 ...
695!> \param nspins ...
696!> \param max_iso_not0 ...
697! **************************************************************************************************
698 SUBROUTINE vh_1c_atom_energy(energy_hartree_1c, ecoul_1c, rho_atom, rrad_0, &
699 grid_atom, iatom, core_charge, rrad_z, Vh1_h, Vh1_s, &
700 nchan_0, nspins, max_iso_not0)
701
702 REAL(dp), INTENT(INOUT) :: energy_hartree_1c
703 TYPE(ecoul_1center_type), DIMENSION(:), POINTER :: ecoul_1c
704 TYPE(rho_atom_type), POINTER :: rho_atom
705 REAL(dp), DIMENSION(:, :), POINTER :: rrad_0
706 TYPE(grid_atom_type), POINTER :: grid_atom
707 INTEGER, INTENT(IN) :: iatom
708 LOGICAL, INTENT(IN) :: core_charge
709 REAL(dp), DIMENSION(:), POINTER :: rrad_z
710 REAL(dp), DIMENSION(:, :), POINTER :: vh1_h, vh1_s
711 INTEGER, INTENT(IN) :: nchan_0, nspins, max_iso_not0
712
713 INTEGER :: iso, ispin, nr
714 REAL(dp) :: ecoul_1_0, ecoul_1_h, ecoul_1_s, &
715 ecoul_1_z
716 TYPE(rho_atom_coeff), DIMENSION(:), POINTER :: r_h, r_s
717
718 nr = grid_atom%nr
719
720 NULLIFY (r_h, r_s)
721 CALL get_rho_atom(rho_atom=rho_atom, rho_rad_h=r_h, rho_rad_s=r_s)
722
723 ! Calculate the contributions to Ecoul coming from Vh1_h*rhoz
724 ecoul_1_z = 0.0_dp
725 IF (core_charge) THEN
726 ecoul_1_z = 0.5_dp*sum(vh1_h(:, 1)*rrad_z(:)*grid_atom%wr(:))
727 END IF
728
729 ! Calculate the contributions to Ecoul coming from Vh1_s*rho0
730 ecoul_1_0 = 0.0_dp
731 DO iso = 1, nchan_0
732 ecoul_1_0 = ecoul_1_0 + 0.5_dp*sum(vh1_s(:, iso)*rrad_0(:, iso)*grid_atom%wr(:))
733 END DO
734
735 ! Calculate the contributions to Ecoul coming from Vh1_h*rho1_h and Vh1_s*rho1_s
736 ecoul_1_s = 0.0_dp
737 ecoul_1_h = 0.0_dp
738 DO ispin = 1, nspins
739 DO iso = 1, max_iso_not0
740 ecoul_1_s = ecoul_1_s + 0.5_dp*sum(vh1_s(:, iso)*r_s(ispin)%r_coef(:, iso)*grid_atom%wr(:))
741 ecoul_1_h = ecoul_1_h + 0.5_dp*sum(vh1_h(:, iso)*r_h(ispin)%r_coef(:, iso)*grid_atom%wr(:))
742 END DO
743 END DO
744
745 CALL set_ecoul_1c(ecoul_1c, iatom, ecoul_1_z=ecoul_1_z, ecoul_1_0=ecoul_1_0)
746 CALL set_ecoul_1c(ecoul_1c=ecoul_1c, iatom=iatom, ecoul_1_h=ecoul_1_h, ecoul_1_s=ecoul_1_s)
747
748 energy_hartree_1c = energy_hartree_1c + ecoul_1_z - ecoul_1_0
749 energy_hartree_1c = energy_hartree_1c + ecoul_1_h - ecoul_1_s
750
751 END SUBROUTINE vh_1c_atom_energy
752
753!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
754
755! **************************************************************************************************
756!> \brief ...
757!> \param rho_atom ...
758!> \param aVh1b_hh ...
759!> \param aVh1b_ss ...
760!> \param aVh1b_00 ...
761!> \param Vh1_h ...
762!> \param Vh1_s ...
763!> \param max_iso_not0 ...
764!> \param max_s_harm ...
765!> \param llmax ...
766!> \param cg_list ...
767!> \param cg_n_list ...
768!> \param nset ...
769!> \param npgf ...
770!> \param lmin ...
771!> \param lmax ...
772!> \param nsotot ...
773!> \param maxso ...
774!> \param nspins ...
775!> \param nchan_0 ...
776!> \param gsph ...
777!> \param g0_h_w ...
778!> \param my_CG ...
779!> \param Qlm_gg ...
780! **************************************************************************************************
781 SUBROUTINE vh_1c_atom_integrals(rho_atom, &
782 aVh1b_hh, aVh1b_ss, aVh1b_00, Vh1_h, Vh1_s, max_iso_not0, &
783 max_s_harm, llmax, cg_list, cg_n_list, &
784 nset, npgf, lmin, lmax, nsotot, maxso, nspins, nchan_0, gsph, &
785 g0_h_w, my_CG, Qlm_gg)
786
787 TYPE(rho_atom_type), POINTER :: rho_atom
788 REAL(dp), DIMENSION(:, :) :: avh1b_hh, avh1b_ss, avh1b_00
789 REAL(dp), DIMENSION(:, :), POINTER :: vh1_h, vh1_s
790 INTEGER, INTENT(IN) :: max_iso_not0, max_s_harm, llmax
791 INTEGER, DIMENSION(:, :, :) :: cg_list
792 INTEGER, DIMENSION(:) :: cg_n_list
793 INTEGER, INTENT(IN) :: nset
794 INTEGER, DIMENSION(:), POINTER :: npgf, lmin, lmax
795 INTEGER, INTENT(IN) :: nsotot, maxso, nspins, nchan_0
796 REAL(dp), DIMENSION(:, :), POINTER :: gsph
797 REAL(dp), DIMENSION(:, 0:) :: g0_h_w
798 REAL(dp), DIMENSION(:, :, :), POINTER :: my_cg, qlm_gg
799
800 INTEGER :: icg, ipgf1, ipgf2, ir, is1, is2, iset1, &
801 iset2, iso, iso1, iso2, ispin, l_ang, &
802 m1, m2, max_iso_not0_local, n1, n2, nr
803 REAL(dp) :: gvg_0, gvg_h, gvg_s
804 TYPE(rho_atom_coeff), DIMENSION(:), POINTER :: int_local_h, int_local_s
805
806 NULLIFY (int_local_h, int_local_s)
807 CALL get_rho_atom(rho_atom=rho_atom, &
808 ga_vlocal_gb_h=int_local_h, &
809 ga_vlocal_gb_s=int_local_s)
810
811 ! Calculate the integrals of the potential with 2 primitives
812 avh1b_hh = 0.0_dp
813 avh1b_ss = 0.0_dp
814 avh1b_00 = 0.0_dp
815
816 nr = SIZE(gsph, 1)
817
818 m1 = 0
819 DO iset1 = 1, nset
820 m2 = 0
821 DO iset2 = 1, nset
822 CALL get_none0_cg_list(my_cg, lmin(iset1), lmax(iset1), lmin(iset2), lmax(iset2), &
823 max_s_harm, llmax, cg_list, cg_n_list, max_iso_not0_local)
824
825 n1 = nsoset(lmax(iset1))
826 DO ipgf1 = 1, npgf(iset1)
827 n2 = nsoset(lmax(iset2))
828 DO ipgf2 = 1, npgf(iset2)
829 ! with contributions to V1_s*rho0
830 DO iso = 1, min(nchan_0, max_iso_not0)
831 l_ang = indso(1, iso)
832 gvg_0 = sum(vh1_s(:, iso)*g0_h_w(:, l_ang))
833 DO icg = 1, cg_n_list(iso)
834 is1 = cg_list(1, icg, iso)
835 is2 = cg_list(2, icg, iso)
836
837 iso1 = is1 + n1*(ipgf1 - 1) + m1
838 iso2 = is2 + n2*(ipgf2 - 1) + m2
839 gvg_h = 0.0_dp
840 gvg_s = 0.0_dp
841
842 DO ir = 1, nr
843 gvg_h = gvg_h + gsph(ir, iso1)*gsph(ir, iso2)*vh1_h(ir, iso)
844 gvg_s = gvg_s + gsph(ir, iso1)*gsph(ir, iso2)*vh1_s(ir, iso)
845 END DO ! ir
846
847 avh1b_hh(iso1, iso2) = avh1b_hh(iso1, iso2) + gvg_h*my_cg(is1, is2, iso)
848 avh1b_ss(iso1, iso2) = avh1b_ss(iso1, iso2) + gvg_s*my_cg(is1, is2, iso)
849 avh1b_00(iso1, iso2) = avh1b_00(iso1, iso2) + gvg_0*qlm_gg(iso1, iso2, iso)
850
851 END DO !icg
852 END DO ! iso
853 ! without contributions to V1_s*rho0
854 DO iso = nchan_0 + 1, max_iso_not0
855 DO icg = 1, cg_n_list(iso)
856 is1 = cg_list(1, icg, iso)
857 is2 = cg_list(2, icg, iso)
858
859 iso1 = is1 + n1*(ipgf1 - 1) + m1
860 iso2 = is2 + n2*(ipgf2 - 1) + m2
861 gvg_h = 0.0_dp
862 gvg_s = 0.0_dp
863
864 DO ir = 1, nr
865 gvg_h = gvg_h + gsph(ir, iso1)*gsph(ir, iso2)*vh1_h(ir, iso)
866 gvg_s = gvg_s + gsph(ir, iso1)*gsph(ir, iso2)*vh1_s(ir, iso)
867 END DO ! ir
868
869 avh1b_hh(iso1, iso2) = avh1b_hh(iso1, iso2) + gvg_h*my_cg(is1, is2, iso)
870 avh1b_ss(iso1, iso2) = avh1b_ss(iso1, iso2) + gvg_s*my_cg(is1, is2, iso)
871
872 END DO !icg
873 END DO ! iso
874 END DO ! ipgf2
875 END DO ! ipgf1
876 m2 = m2 + maxso
877 END DO ! iset2
878 m1 = m1 + maxso
879 END DO !iset1
880 DO ispin = 1, nspins
881 CALL daxpy(nsotot*nsotot, 1.0_dp, avh1b_hh, 1, int_local_h(ispin)%r_coef, 1)
882 CALL daxpy(nsotot*nsotot, 1.0_dp, avh1b_ss, 1, int_local_s(ispin)%r_coef, 1)
883 CALL daxpy(nsotot*nsotot, -1.0_dp, avh1b_00, 1, int_local_h(ispin)%r_coef, 1)
884 CALL daxpy(nsotot*nsotot, -1.0_dp, avh1b_00, 1, int_local_s(ispin)%r_coef, 1)
885 END DO ! ispin
886
887 END SUBROUTINE vh_1c_atom_integrals
888
889!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
890
891END MODULE hartree_local_methods
Define the atomic kind types and their sub types.
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.
subroutine, public get_gto_basis_set(gto_basis_set, name, aliases, norm_type, kind_radius, ncgf, nset, nsgf, cgf_symbol, sgf_symbol, norm_cgf, set_radius, lmax, lmin, lx, ly, lz, m, ncgf_set, npgf, nsgf_set, nshell, cphi, pgf_radius, sphi, scon, zet, first_cgf, first_sgf, l, last_cgf, last_sgf, n, gcc, maxco, maxl, maxpgf, maxsgf_set, maxshell, maxso, nco_sum, npgf_sum, nshell_sum, maxder, short_kind_radius, npgf_seg_sum, ccon)
...
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 init_coulomb_local(hartree_local, natom)
...
subroutine, public vh_1c_gg_integrals(qs_env, energy_hartree_1c, ecoul_1c, local_rho_set, para_env, tddft, local_rho_set_2nd, core_2nd)
Calculates one center GAPW Hartree energies and matrix elements Hartree potentials are input Takes po...
subroutine, public calculate_vh_1center(vrad_h, vrad_s, rrad_h, rrad_s, rrad_0, rrad_z, grid_atom)
Calculates Hartree potential for hard and soft densities (including nuclear charge and compensation c...
subroutine, public allocate_ecoul_1center(ecoul_1c, natom)
...
subroutine, public set_ecoul_1c(ecoul_1c, iatom, ecoul_1_h, ecoul_1_s, ecoul_1_z, ecoul_1_0)
...
Defines the basic variable types.
Definition kinds.F:23
integer, parameter, public dp
Definition kinds.F:34
Definition of mathematical constants and functions.
real(kind=dp), parameter, public pi
real(kind=dp), parameter, public fourpi
Interface to the message passing library MPI.
Provides Cartesian and spherical orbital pointers and indices.
integer, dimension(:), allocatable, public nsoset
integer, dimension(:, :), allocatable, public indso
container for various plainwaves related things
subroutine, public pw_env_get(pw_env, pw_pools, cube_info, gridlevel_info, auxbas_pw_pool, auxbas_grid, auxbas_rs_desc, auxbas_rs_grid, rs_descs, rs_grids, xc_pw_pool, vdw_pw_pool, poisson_env, interp_section)
returns the various attributes of the pw env
functions related to the poisson solver on regular grids
integer, parameter, public pw_poisson_periodic
container for information about total charges on the grids
A collection of functions used by CNEO-DFT (see J. Chem. Theory Comput. 2025, 21, 16,...
subroutine, public vh_1c_nuc_integrals(rhoz_cneo, zeff, avh1b_hh, avh1b_ss, avh1b_00, vh1_h, vh1_s, max_iso_not0_elec, max_iso_not0_nuc, max_s_harm, llmax, cg_list, cg_n_list, nset, npgf, lmin, lmax, nsotot, maxso, nchan_0, gsph, g0_h_w, my_cg, qlm_gg)
Mostly copied from hartree_local_methods::Vh_1c_atom_integrals.
subroutine, public calculate_rhoz_cneo(rho, potential, cg_list, cg_n_list, nset, npgf, lmin, lmax, maxl, maxso)
...
Types used by CNEO-DFT (see J. Chem. Theory Comput. 2025, 21, 16, 7865–7877)
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, matrix_vhxc, 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.
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, hund_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, proj_shell_charge, lr_atom, do_mtlr, u_j_loop, ao_coef, 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.
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, do_mtlr_present, basis_type, total_zeff_corr, npgf_seg, cneo_potential_present, nkind_q, natom_q)
Get attributes of an atomic kind set.
subroutine, public get_local_rho(local_rho_set, rho_atom_set, rho0_atom_set, rho0_mpole, rhoz_set, rhoz_cneo_set)
...
subroutine, public get_rho0_mpole(rho0_mpole, g0_h, vg0_h, iat, ikind, lmax_0, l0_ikind, mp_gau_ikind, mp_rho, norm_g0l_h, qlm_gg, qlm_car, qlm_tot, zet0_h, igrid_zet0_s, rpgf0_h, rpgf0_s, max_rpgf0_s, rho0_s_rs, rho0_s_gs, rhoz_cneo_s_rs, rhoz_cneo_s_gs)
...
subroutine, public get_rho_atom(rho_atom, cpc_h, cpc_s, rho_rad_h, rho_rad_s, drho_rad_h, drho_rad_s, vrho_rad_h, vrho_rad_s, rho_rad_h_d, rho_rad_s_d, ga_vlocal_gb_h, ga_vlocal_gb_s, int_scr_h, int_scr_s)
...
All kind of helpful little routines.
Definition util.F:14
pure integer function, dimension(2), public get_limit(m, n, me)
divide m entries into n parts, return size of part me
Definition util.F:333
Provides all information about an atomic kind.
Type defining parameters related to the simulation cell.
Definition cell_types.F:60
stores all the informations relevant to an mpi environment
contained for different pw related things
environment for the poisson solver
Container for information about total charges on the grids.
Provides all information about a quickstep kind.