(git:b17b328)
Loading...
Searching...
No Matches
qs_rho0_types.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! **************************************************************************************************
8MODULE qs_rho0_types
9
10 USE cp_units, ONLY: cp_unit_from_cp2k
11 USE kinds, ONLY: dp
12 USE mathconstants, ONLY: fourpi,&
13 pi,&
14 rootpi
15 USE memory_utilities, ONLY: reallocate
16 USE pw_types, ONLY: pw_c1d_gs_type,&
17 pw_r3d_rs_type
18 USE qs_grid_atom, ONLY: grid_atom_type
19 USE qs_rho_atom_types, ONLY: rho_atom_coeff
20 USE whittaker, ONLY: whittaker_c0a,&
21 whittaker_ci
22#include "./base/base_uses.f90"
23
24 IMPLICIT NONE
25
26 PRIVATE
27
28! *** Global parameters (only in this module)
29
30 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'qs_rho0_types'
31
32! *** Define multipole type ***
33
34! **************************************************************************************************
35 TYPE mpole_rho_atom
36 REAL(dp), DIMENSION(:), POINTER :: qlm_h => null(), &
37 qlm_s => null(), &
38 qlm_tot => null(), &
39 qlm_car => null()
40 REAL(dp) :: qlm_z = -1.0_dp
41 REAL(dp), DIMENSION(2) :: q0 = -1.0_dp
42 END TYPE mpole_rho_atom
43
44! **************************************************************************************************
45 TYPE mpole_gau_overlap
46 REAL(dp), DIMENSION(:, :, :), POINTER :: qlm_gg => null()
47 REAL(dp), DIMENSION(:, :), POINTER :: g0_h => null(), vg0_h => null()
48 REAL(dp) :: rpgf0_h = -1.0_dp, rpgf0_s = -1.0_dp
49 END TYPE mpole_gau_overlap
50
51! **************************************************************************************************
52 TYPE rho0_mpole_type
53 TYPE(mpole_rho_atom), DIMENSION(:), POINTER :: mp_rho => null()
54 TYPE(mpole_gau_overlap), DIMENSION(:), &
55 POINTER :: mp_gau => null()
56 REAL(dp) :: zet0_h = -1.0_dp, &
57 total_rho0_h = -1.0_dp
58 REAL(dp) :: max_rpgf0_s = -1.0_dp
59 REAL(dp), DIMENSION(:), POINTER :: norm_g0l_h => null()
60 INTEGER, DIMENSION(:), POINTER :: lmax0_kind => null()
61 INTEGER :: lmax_0 = -1, igrid_zet0_s = -1
62 TYPE(pw_r3d_rs_type), POINTER :: rho0_s_rs => null()
63 TYPE(pw_c1d_gs_type), POINTER :: rho0_s_gs => null()
64
65 ! CNEO nuclear charge density related stuff.
66 ! These are put here because it is preferred that function rho0_s_grid_create
67 ! can initialize PW for both rho0 and nuclear rho1s. This is akin to the behavior
68 ! that init_rho0 also initializes data structures for CNEO nuclear charge densities,
69 ! such that code modifications are mimimal for places where init_rho0 and
70 ! rho0_s_grid_create are called.
71 LOGICAL :: do_cneo = .false. ! if CNEO potential is used
72 REAL(dp) :: tot_rhoz_cneo_s = 0.0_dp ! soft nuclear charge on the local grid
73 TYPE(pw_r3d_rs_type), POINTER :: rhoz_cneo_s_rs => null() ! soft nuclear charge density, r-space
74 TYPE(pw_c1d_gs_type), POINTER :: rhoz_cneo_s_gs => null() ! soft nuclear charge density, g-space
75 END TYPE rho0_mpole_type
76
77! **************************************************************************************************
78 TYPE rho0_atom_type
79 TYPE(rho_atom_coeff), POINTER :: rho0_rad_h => null(), &
80 vrho0_rad_h => null()
81 END TYPE rho0_atom_type
82
83! Public Types
84
85 PUBLIC :: mpole_rho_atom, mpole_gau_overlap, &
86 rho0_atom_type, rho0_mpole_type
87
88! Public Subroutine
89
90 PUBLIC :: allocate_multipoles, allocate_rho0_mpole, &
91 allocate_rho0_atom, allocate_rho0_atom_rad, &
92 deallocate_rho0_atom, deallocate_rho0_mpole, &
93 calculate_g0, get_rho0_mpole, initialize_mpole_rho, &
94 write_rho0_info
95
96CONTAINS
97
98! **************************************************************************************************
99!> \brief ...
100!> \param mp_rho ...
101!> \param natom ...
102!> \param mp_gau ...
103!> \param nkind ...
104! **************************************************************************************************
105 SUBROUTINE allocate_multipoles(mp_rho, natom, mp_gau, nkind)
106
107 TYPE(mpole_rho_atom), DIMENSION(:), POINTER :: mp_rho
108 INTEGER, INTENT(IN) :: natom
109 TYPE(mpole_gau_overlap), DIMENSION(:), POINTER :: mp_gau
110 INTEGER, INTENT(IN) :: nkind
111
112 INTEGER :: iat, ikind
113
114 IF (ASSOCIATED(mp_rho)) THEN
115 CALL deallocate_mpole_rho(mp_rho)
116 END IF
117
118 ALLOCATE (mp_rho(natom))
119
120 DO iat = 1, natom
121 NULLIFY (mp_rho(iat)%Qlm_h)
122 NULLIFY (mp_rho(iat)%Qlm_s)
123 NULLIFY (mp_rho(iat)%Qlm_tot)
124 NULLIFY (mp_rho(iat)%Qlm_car)
125 END DO
126
127 IF (ASSOCIATED(mp_gau)) THEN
128 CALL deallocate_mpole_gau(mp_gau)
129 END IF
130
131 ALLOCATE (mp_gau(nkind))
132
133 DO ikind = 1, nkind
134 NULLIFY (mp_gau(ikind)%Qlm_gg)
135 NULLIFY (mp_gau(ikind)%g0_h)
136 NULLIFY (mp_gau(ikind)%vg0_h)
137 mp_gau(ikind)%rpgf0_h = 0.0_dp
138 mp_gau(ikind)%rpgf0_s = 0.0_dp
139 END DO
140
141 END SUBROUTINE allocate_multipoles
142
143! **************************************************************************************************
144!> \brief ...
145!> \param rho0_set ...
146!> \param natom ...
147! **************************************************************************************************
148 SUBROUTINE allocate_rho0_atom(rho0_set, natom)
149
150 TYPE(rho0_atom_type), DIMENSION(:), POINTER :: rho0_set
151 INTEGER, INTENT(IN) :: natom
152
153 INTEGER :: iat
154
155 IF (ASSOCIATED(rho0_set)) THEN
156 CALL deallocate_rho0_atom(rho0_set)
157 END IF
158
159 ALLOCATE (rho0_set(natom))
160
161 DO iat = 1, natom
162 NULLIFY (rho0_set(iat)%rho0_rad_h)
163 NULLIFY (rho0_set(iat)%vrho0_rad_h)
164 END DO
165
166 END SUBROUTINE allocate_rho0_atom
167
168! **************************************************************************************************
169!> \brief ...
170!> \param rho0_atom ...
171!> \param nr ...
172!> \param nchannels ...
173! **************************************************************************************************
174 SUBROUTINE allocate_rho0_atom_rad(rho0_atom, nr, nchannels)
175
176 TYPE(rho0_atom_type), INTENT(OUT) :: rho0_atom
177 INTEGER, INTENT(IN) :: nr, nchannels
178
179 ALLOCATE (rho0_atom%rho0_rad_h)
180
181 NULLIFY (rho0_atom%rho0_rad_h%r_coef)
182 ALLOCATE (rho0_atom%rho0_rad_h%r_coef(1:nr, 1:nchannels))
183 rho0_atom%rho0_rad_h%r_coef = 0.0_dp
184
185 ALLOCATE (rho0_atom%vrho0_rad_h)
186
187 NULLIFY (rho0_atom%vrho0_rad_h%r_coef)
188 ALLOCATE (rho0_atom%vrho0_rad_h%r_coef(1:nr, 1:nchannels))
189 rho0_atom%vrho0_rad_h%r_coef = 0.0_dp
190
191 END SUBROUTINE allocate_rho0_atom_rad
192
193! **************************************************************************************************
194!> \brief ...
195!> \param rho0 ...
196! **************************************************************************************************
197 SUBROUTINE allocate_rho0_mpole(rho0)
198
199 TYPE(rho0_mpole_type), POINTER :: rho0
200
201 IF (ASSOCIATED(rho0)) THEN
202 CALL deallocate_rho0_mpole(rho0)
203 END IF
204
205 ALLOCATE (rho0)
206
207 NULLIFY (rho0%mp_rho)
208 NULLIFY (rho0%mp_gau)
209 NULLIFY (rho0%norm_g0l_h)
210 NULLIFY (rho0%lmax0_kind)
211 NULLIFY (rho0%rho0_s_rs)
212 NULLIFY (rho0%rho0_s_gs)
213
214 END SUBROUTINE allocate_rho0_mpole
215
216! **************************************************************************************************
217!> \brief ...
218!> \param rho0_mpole ...
219!> \param grid_atom ...
220!> \param ik ...
221! **************************************************************************************************
222 SUBROUTINE calculate_g0(rho0_mpole, grid_atom, ik)
223
224 TYPE(rho0_mpole_type), POINTER :: rho0_mpole
225 TYPE(grid_atom_type), POINTER :: grid_atom
226 INTEGER, INTENT(IN) :: ik
227
228 INTEGER :: ir, l, lmax, nr
229 REAL(dp) :: c1, prefactor, root_z_h, z_h
230 REAL(dp), ALLOCATABLE, DIMENSION(:) :: erf_z_h, gexp, gh_tmp, int1, int2
231
232 nr = grid_atom%nr
233 lmax = rho0_mpole%lmax0_kind(ik)
234 z_h = rho0_mpole%zet0_h
235 root_z_h = sqrt(z_h)
236
237! Allocate g0
238 CALL reallocate(rho0_mpole%mp_gau(ik)%g0_h, 1, nr, 0, lmax)
239 CALL reallocate(rho0_mpole%mp_gau(ik)%vg0_h, 1, nr, 0, lmax)
240
241 ALLOCATE (gexp(nr), gh_tmp(nr), erf_z_h(nr), int1(nr), int2(nr))
242
243 gh_tmp(1:nr) = exp(-z_h*grid_atom%rad2(1:nr))
244
245 DO ir = 1, nr
246 erf_z_h(ir) = erf(grid_atom%rad(ir)*root_z_h)
247 END DO
248
249 DO ir = 1, nr
250 IF (gh_tmp(ir) < 1.0e-30_dp) gh_tmp(ir) = 0.0_dp
251 END DO
252
253 gexp(1:nr) = gh_tmp(1:nr)
254 rho0_mpole%mp_gau(ik)%g0_h(1:nr, 0) = gh_tmp(1:nr)* &
255 rho0_mpole%norm_g0l_h(0)
256 CALL whittaker_c0a(int1, grid_atom%rad, gh_tmp, erf_z_h, z_h, 0, 0, nr)
257 CALL whittaker_ci(int2, grid_atom%rad, gh_tmp, z_h, 0, nr)
258
259 prefactor = fourpi*rho0_mpole%norm_g0l_h(0)
260
261 c1 = sqrt(pi*pi*pi/(z_h*z_h*z_h))*rho0_mpole%norm_g0l_h(0)
262
263 DO ir = 1, nr
264 rho0_mpole%mp_gau(ik)%vg0_h(ir, 0) = c1*erf_z_h(ir)*grid_atom%oorad2l(ir, 1)
265 END DO
266
267 DO l = 1, lmax
268 gh_tmp(1:nr) = gh_tmp(1:nr)*grid_atom%rad(1:nr)
269 rho0_mpole%mp_gau(ik)%g0_h(1:nr, l) = gh_tmp(1:nr)* &
270 rho0_mpole%norm_g0l_h(l)
271
272 prefactor = fourpi/(2.0_dp*l + 1.0_dp)*rho0_mpole%norm_g0l_h(l)
273 CALL whittaker_c0a(int1, grid_atom%rad, gexp, erf_z_h, z_h, l, l, nr)
274 DO ir = 1, nr
275 rho0_mpole%mp_gau(ik)%vg0_h(ir, l) = prefactor*(int1(ir) + &
276 int2(ir)*grid_atom%rad2l(ir, l))
277 END DO
278
279 END DO ! l
280
281 DEALLOCATE (gexp, erf_z_h, gh_tmp, int1, int2)
282 END SUBROUTINE calculate_g0
283
284! **************************************************************************************************
285!> \brief ...
286!> \param mp_gau ...
287! **************************************************************************************************
288 SUBROUTINE deallocate_mpole_gau(mp_gau)
289
290 TYPE(mpole_gau_overlap), DIMENSION(:), POINTER :: mp_gau
291
292 INTEGER :: ikind, nkind
293
294 nkind = SIZE(mp_gau)
295
296 DO ikind = 1, nkind
297
298 IF (ASSOCIATED(mp_gau(ikind)%Qlm_gg)) THEN
299 DEALLOCATE (mp_gau(ikind)%Qlm_gg)
300 END IF
301
302 DEALLOCATE (mp_gau(ikind)%g0_h)
303
304 DEALLOCATE (mp_gau(ikind)%vg0_h)
305 END DO
306
307 DEALLOCATE (mp_gau)
308
309 END SUBROUTINE deallocate_mpole_gau
310
311! **************************************************************************************************
312!> \brief ...
313!> \param mp_rho ...
314! **************************************************************************************************
315 SUBROUTINE deallocate_mpole_rho(mp_rho)
316
317 TYPE(mpole_rho_atom), DIMENSION(:), POINTER :: mp_rho
318
319 INTEGER :: iat, natom
320
321 natom = SIZE(mp_rho)
322
323 DO iat = 1, natom
324 DEALLOCATE (mp_rho(iat)%Qlm_h)
325 DEALLOCATE (mp_rho(iat)%Qlm_s)
326 DEALLOCATE (mp_rho(iat)%Qlm_tot)
327 DEALLOCATE (mp_rho(iat)%Qlm_car)
328 END DO
329
330 DEALLOCATE (mp_rho)
331
332 END SUBROUTINE deallocate_mpole_rho
333
334! **************************************************************************************************
335!> \brief ...
336!> \param rho0_atom_set ...
337! **************************************************************************************************
338 SUBROUTINE deallocate_rho0_atom(rho0_atom_set)
339
340 TYPE(rho0_atom_type), DIMENSION(:), POINTER :: rho0_atom_set
341
342 INTEGER :: iat, natom
343
344 IF (ASSOCIATED(rho0_atom_set)) THEN
345
346 natom = SIZE(rho0_atom_set)
347
348 DO iat = 1, natom
349 IF (ASSOCIATED(rho0_atom_set(iat)%rho0_rad_h)) THEN
350 DEALLOCATE (rho0_atom_set(iat)%rho0_rad_h%r_coef)
351 DEALLOCATE (rho0_atom_set(iat)%rho0_rad_h)
352 END IF
353 IF (ASSOCIATED(rho0_atom_set(iat)%vrho0_rad_h)) THEN
354 DEALLOCATE (rho0_atom_set(iat)%vrho0_rad_h%r_coef)
355 DEALLOCATE (rho0_atom_set(iat)%vrho0_rad_h)
356 END IF
357 END DO
358
359 DEALLOCATE (rho0_atom_set)
360 ELSE
361 CALL cp_abort(__location__, &
362 "The pointer rho0_atom_set is not associated and "// &
363 "cannot be deallocated")
364 END IF
365
366 END SUBROUTINE deallocate_rho0_atom
367! **************************************************************************************************
368!> \brief ...
369!> \param rho0 ...
370! **************************************************************************************************
371 SUBROUTINE deallocate_rho0_mpole(rho0)
372
373 TYPE(rho0_mpole_type), POINTER :: rho0
374
375 IF (ASSOCIATED(rho0)) THEN
376
377 IF (ASSOCIATED(rho0%mp_gau)) CALL deallocate_mpole_gau(rho0%mp_gau)
378
379 IF (ASSOCIATED(rho0%mp_rho)) CALL deallocate_mpole_rho(rho0%mp_rho)
380
381 IF (ASSOCIATED(rho0%lmax0_kind)) THEN
382 DEALLOCATE (rho0%lmax0_kind)
383 END IF
384
385 IF (ASSOCIATED(rho0%norm_g0l_h)) THEN
386 DEALLOCATE (rho0%norm_g0l_h)
387 END IF
388
389 IF (ASSOCIATED(rho0%rho0_s_rs)) THEN
390 CALL rho0%rho0_s_rs%release()
391 DEALLOCATE (rho0%rho0_s_rs)
392 END IF
393
394 IF (ASSOCIATED(rho0%rho0_s_gs)) THEN
395 CALL rho0%rho0_s_gs%release()
396 DEALLOCATE (rho0%rho0_s_gs)
397 END IF
398
399 IF (ASSOCIATED(rho0%rhoz_cneo_s_rs)) THEN
400 CALL rho0%rhoz_cneo_s_rs%release()
401 DEALLOCATE (rho0%rhoz_cneo_s_rs)
402 END IF
403
404 IF (ASSOCIATED(rho0%rhoz_cneo_s_gs)) THEN
405 CALL rho0%rhoz_cneo_s_gs%release()
406 DEALLOCATE (rho0%rhoz_cneo_s_gs)
407 END IF
408
409 DEALLOCATE (rho0)
410 ELSE
411 CALL cp_abort(__location__, &
412 "The pointer rho0 is not associated and "// &
413 "cannot be deallocated")
414 END IF
415
416 END SUBROUTINE deallocate_rho0_mpole
417
418! **************************************************************************************************
419!> \brief ...
420!> \param rho0_mpole ...
421!> \param g0_h ...
422!> \param vg0_h ...
423!> \param iat ...
424!> \param ikind ...
425!> \param lmax_0 ...
426!> \param l0_ikind ...
427!> \param mp_gau_ikind ...
428!> \param mp_rho ...
429!> \param norm_g0l_h ...
430!> \param Qlm_gg ...
431!> \param Qlm_car ...
432!> \param Qlm_tot ...
433!> \param zet0_h ...
434!> \param igrid_zet0_s ...
435!> \param rpgf0_h ...
436!> \param rpgf0_s ...
437!> \param max_rpgf0_s ...
438!> \param rho0_s_rs ...
439!> \param rho0_s_gs ...
440!> \param rhoz_cneo_s_rs ...
441!> \param rhoz_cneo_s_gs ...
442! **************************************************************************************************
443 SUBROUTINE get_rho0_mpole(rho0_mpole, g0_h, vg0_h, iat, ikind, lmax_0, l0_ikind, &
444 mp_gau_ikind, mp_rho, norm_g0l_h, &
445 Qlm_gg, Qlm_car, Qlm_tot, &
446 zet0_h, igrid_zet0_s, rpgf0_h, rpgf0_s, &
447 max_rpgf0_s, rho0_s_rs, rho0_s_gs, &
448 rhoz_cneo_s_rs, rhoz_cneo_s_gs)
449
450 TYPE(rho0_mpole_type), POINTER :: rho0_mpole
451 REAL(dp), DIMENSION(:, :), OPTIONAL, POINTER :: g0_h, vg0_h
452 INTEGER, INTENT(IN), OPTIONAL :: iat, ikind
453 INTEGER, INTENT(OUT), OPTIONAL :: lmax_0, l0_ikind
454 TYPE(mpole_gau_overlap), OPTIONAL, POINTER :: mp_gau_ikind
455 TYPE(mpole_rho_atom), DIMENSION(:), OPTIONAL, &
456 POINTER :: mp_rho
457 REAL(dp), DIMENSION(:), OPTIONAL, POINTER :: norm_g0l_h
458 REAL(dp), DIMENSION(:, :, :), OPTIONAL, POINTER :: qlm_gg
459 REAL(dp), DIMENSION(:), OPTIONAL, POINTER :: qlm_car, qlm_tot
460 REAL(dp), INTENT(OUT), OPTIONAL :: zet0_h
461 INTEGER, INTENT(OUT), OPTIONAL :: igrid_zet0_s
462 REAL(dp), INTENT(OUT), OPTIONAL :: rpgf0_h, rpgf0_s, max_rpgf0_s
463 TYPE(pw_r3d_rs_type), OPTIONAL, POINTER :: rho0_s_rs
464 TYPE(pw_c1d_gs_type), OPTIONAL, POINTER :: rho0_s_gs
465 TYPE(pw_r3d_rs_type), OPTIONAL, POINTER :: rhoz_cneo_s_rs
466 TYPE(pw_c1d_gs_type), OPTIONAL, POINTER :: rhoz_cneo_s_gs
467
468 IF (ASSOCIATED(rho0_mpole)) THEN
469
470 IF (PRESENT(lmax_0)) lmax_0 = rho0_mpole%lmax_0
471 IF (PRESENT(mp_rho)) mp_rho => rho0_mpole%mp_rho
472 IF (PRESENT(norm_g0l_h)) norm_g0l_h => rho0_mpole%norm_g0l_h
473 IF (PRESENT(zet0_h)) zet0_h = rho0_mpole%zet0_h
474 IF (PRESENT(igrid_zet0_s)) igrid_zet0_s = rho0_mpole%igrid_zet0_s
475 IF (PRESENT(max_rpgf0_s)) max_rpgf0_s = rho0_mpole%max_rpgf0_s
476 IF (PRESENT(rho0_s_rs)) rho0_s_rs => rho0_mpole%rho0_s_rs
477 IF (PRESENT(rho0_s_gs)) rho0_s_gs => rho0_mpole%rho0_s_gs
478 IF (PRESENT(rhoz_cneo_s_rs)) rhoz_cneo_s_rs => rho0_mpole%rhoz_cneo_s_rs
479 IF (PRESENT(rhoz_cneo_s_gs)) rhoz_cneo_s_gs => rho0_mpole%rhoz_cneo_s_gs
480
481 IF (PRESENT(ikind)) THEN
482 IF (PRESENT(l0_ikind)) l0_ikind = rho0_mpole%lmax0_kind(ikind)
483 IF (PRESENT(mp_gau_ikind)) mp_gau_ikind => rho0_mpole%mp_gau(ikind)
484 IF (PRESENT(g0_h)) g0_h => rho0_mpole%mp_gau(ikind)%g0_h
485 IF (PRESENT(vg0_h)) vg0_h => rho0_mpole%mp_gau(ikind)%vg0_h
486 IF (PRESENT(qlm_gg)) qlm_gg => rho0_mpole%mp_gau(ikind)%Qlm_gg
487 IF (PRESENT(rpgf0_h)) rpgf0_h = rho0_mpole%mp_gau(ikind)%rpgf0_h
488 IF (PRESENT(rpgf0_s)) rpgf0_s = rho0_mpole%mp_gau(ikind)%rpgf0_s
489 END IF
490 IF (PRESENT(iat)) THEN
491 IF (PRESENT(qlm_car)) qlm_car => rho0_mpole%mp_rho(iat)%Qlm_car
492 IF (PRESENT(qlm_tot)) qlm_tot => rho0_mpole%mp_rho(iat)%Qlm_tot
493 END IF
494
495 ELSE
496 cpabort("The pointer rho0_mpole is not associated")
497 END IF
498
499 END SUBROUTINE get_rho0_mpole
500
501! **************************************************************************************************
502!> \brief ...
503!> \param mp_rho ...
504!> \param nchan_s ...
505!> \param nchan_c ...
506!> \param zeff ...
507! **************************************************************************************************
508 SUBROUTINE initialize_mpole_rho(mp_rho, nchan_s, nchan_c, zeff)
509
510 TYPE(mpole_rho_atom) :: mp_rho
511 INTEGER, INTENT(IN) :: nchan_s, nchan_c
512 REAL(kind=dp), INTENT(IN) :: zeff
513
514 CALL reallocate(mp_rho%Qlm_h, 1, nchan_s)
515 CALL reallocate(mp_rho%Qlm_s, 1, nchan_s)
516 CALL reallocate(mp_rho%Qlm_tot, 1, nchan_s)
517 CALL reallocate(mp_rho%Qlm_car, 1, nchan_c)
518
519 mp_rho%Qlm_h = 0.0_dp
520 mp_rho%Qlm_s = 0.0_dp
521 mp_rho%Qlm_tot = 0.0_dp
522 mp_rho%Qlm_car = 0.0_dp
523 mp_rho%Qlm_z = -2.0_dp*rootpi*zeff
524 mp_rho%Q0 = 0.0_dp
525
526 END SUBROUTINE initialize_mpole_rho
527
528! **************************************************************************************************
529!> \brief ...
530!> \param rho0_mpole ...
531!> \param unit_str ...
532!> \param output_unit ...
533! **************************************************************************************************
534 SUBROUTINE write_rho0_info(rho0_mpole, unit_str, output_unit)
535
536 TYPE(rho0_mpole_type), POINTER :: rho0_mpole
537 CHARACTER(LEN=*), INTENT(IN) :: unit_str
538 INTEGER, INTENT(in) :: output_unit
539
540 INTEGER :: ikind, l, nkind
541 REAL(dp) :: conv
542
543 IF (ASSOCIATED(rho0_mpole)) THEN
544 conv = cp_unit_from_cp2k(1.0_dp, trim(unit_str))
545
546 WRITE (unit=output_unit, fmt="(/,T5,A,/)") &
547 "*** Compensation density charges data set ***"
548 WRITE (unit=output_unit, fmt="(T2,A,T35,f16.10)") &
549 "- Rho0 exponent :", rho0_mpole%zet0_h
550 WRITE (unit=output_unit, fmt="(T2,A,T35,I5)") &
551 "- Global max l :", rho0_mpole%lmax_0
552
553 WRITE (unit=output_unit, fmt="(T2,A)") &
554 "- Normalization constants for g0"
555 DO l = 0, rho0_mpole%lmax_0
556 WRITE (unit=output_unit, fmt="(T20,A,T31,I2,T38,A,f15.5)") &
557 "ang. mom.= ", l, " hard= ", rho0_mpole%norm_g0l_h(l)
558 END DO
559
560 nkind = SIZE(rho0_mpole%lmax0_kind, 1)
561 DO ikind = 1, nkind
562 WRITE (unit=output_unit, fmt="(/,T2,A,T55,I2)") &
563 "- rho0 max L and radii in "//trim(unit_str)// &
564 " for the atom kind :", ikind
565
566 WRITE (unit=output_unit, fmt="(T2,T20,A,T55,I5)") &
567 "=> l max :", rho0_mpole%lmax0_kind(ikind)
568
569 WRITE (unit=output_unit, fmt="(T2,T20,A,T55,f20.10)") &
570 "=> max radius of g0: ", &
571 rho0_mpole%mp_gau(ikind)%rpgf0_h*conv
572 END DO ! ikind
573
574 ELSE
575 WRITE (unit=output_unit, fmt="(/,T5,A,/)") &
576 ' WARNING: I cannot print rho0, it is not associated'
577 END IF
578
579 END SUBROUTINE write_rho0_info
580END MODULE qs_rho0_types
memory_utilities::reallocate
Definition memory_utilities.F:27
cp_units
unit conversion facility
Definition cp_units.F:30
cp_units::cp_unit_from_cp2k
real(kind=dp) function, public cp_unit_from_cp2k(value, unit_str, defaults, power)
converts from the internal cp2k units to the given unit
Definition cp_units.F:1178
kinds
Defines the basic variable types.
Definition kinds.F:23
kinds::dp
integer, parameter, public dp
Definition kinds.F:34
mathconstants
Definition of mathematical constants and functions.
Definition mathconstants.F:16
mathconstants::pi
real(kind=dp), parameter, public pi
Definition mathconstants.F:110
mathconstants::rootpi
real(kind=dp), parameter, public rootpi
Definition mathconstants.F:114
mathconstants::fourpi
real(kind=dp), parameter, public fourpi
Definition mathconstants.F:113
memory_utilities
Utility routines for the memory handling.
Definition memory_utilities.F:14
pw_types
Definition pw_types.F:24
qs_grid_atom
Definition qs_grid_atom.F:8
qs_rho0_types
Definition qs_rho0_types.F:8
qs_rho0_types::allocate_multipoles
subroutine, public allocate_multipoles(mp_rho, natom, mp_gau, nkind)
...
Definition qs_rho0_types.F:106
qs_rho0_types::initialize_mpole_rho
subroutine, public initialize_mpole_rho(mp_rho, nchan_s, nchan_c, zeff)
...
Definition qs_rho0_types.F:509
qs_rho0_types::allocate_rho0_atom_rad
subroutine, public allocate_rho0_atom_rad(rho0_atom, nr, nchannels)
...
Definition qs_rho0_types.F:175
qs_rho0_types::calculate_g0
subroutine, public calculate_g0(rho0_mpole, grid_atom, ik)
...
Definition qs_rho0_types.F:223
qs_rho0_types::write_rho0_info
subroutine, public write_rho0_info(rho0_mpole, unit_str, output_unit)
...
Definition qs_rho0_types.F:535
qs_rho0_types::allocate_rho0_atom
subroutine, public allocate_rho0_atom(rho0_set, natom)
...
Definition qs_rho0_types.F:149
qs_rho0_types::allocate_rho0_mpole
subroutine, public allocate_rho0_mpole(rho0)
...
Definition qs_rho0_types.F:198
qs_rho0_types::deallocate_rho0_mpole
subroutine, public deallocate_rho0_mpole(rho0)
...
Definition qs_rho0_types.F:372
qs_rho0_types::deallocate_rho0_atom
subroutine, public deallocate_rho0_atom(rho0_atom_set)
...
Definition qs_rho0_types.F:339
qs_rho0_types::get_rho0_mpole
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)
...
Definition qs_rho0_types.F:449
qs_rho_atom_types
Definition qs_rho_atom_types.F:8
whittaker
Calculates special integrals.
Definition whittaker.F:12
whittaker::whittaker_c0a
subroutine, public whittaker_c0a(wc, r, expa, erfa, alpha, l1, l2, n)
int(y^(2+l1+l2) * exp(-alpha*y*y),y=0..x) / x^(l2+1); wc(:) :: output r(:) :: coordinate expa(:) :: e...
Definition whittaker.F:52
whittaker::whittaker_ci
subroutine, public whittaker_ci(wc, r, expa, alpha, l, n)
int(y^(l+1) * exp(-alpha*y*y),y=x..infinity);
Definition whittaker.F:340
pw_types::pw_c1d_gs_type
Definition pw_types.F:127
pw_types::pw_r3d_rs_type
Definition pw_types.F:62
qs_grid_atom::grid_atom_type
Definition qs_grid_atom.F:44
qs_rho0_types::mpole_gau_overlap
Definition qs_rho0_types.F:45
qs_rho0_types::mpole_rho_atom
Definition qs_rho0_types.F:35
qs_rho0_types::rho0_atom_type
Definition qs_rho0_types.F:78
qs_rho0_types::rho0_mpole_type
Definition qs_rho0_types.F:52
qs_rho_atom_types::rho_atom_coeff
Definition qs_rho_atom_types.F:19