(git:81b4ad0)
Loading...
Searching...
No Matches
qs_interactions.F
Go to the documentation of this file.
1!--------------------------------------------------------------------------------------------------!
2! CP2K: A general program to perform molecular dynamics simulations !
3! Copyright 2000-2026 CP2K developers group <https://cp2k.org> !
4! !
5! SPDX-License-Identifier: GPL-2.0-or-later !
6!--------------------------------------------------------------------------------------------------!
7
8! **************************************************************************************************
9!> \brief Calculate the interaction radii for the operator matrix calculation.
10!> \par History
11!> Joost VandeVondele : added exp_radius_very_extended
12!> 24.09.2002 overloaded init_interaction_radii for KG use (gt)
13!> 27.02.2004 integrated KG into QS version of init_interaction_radii (jgh)
14!> 07.03.2006 new routine to extend kind interaction radius for semi-empiric (jgh)
15!> \author MK (12.07.2000)
16! **************************************************************************************************
18
19 USE ao_util, ONLY: exp_radius
39 USE kinds, ONLY: default_string_length,&
40 dp
45 USE qs_kind_types, ONLY: get_qs_kind,&
49#include "./base/base_uses.f90"
50
51 IMPLICIT NONE
52
53 PRIVATE
54
55! *** Global parameters (only in this module)
56
57 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'qs_interactions'
58
59! *** Public subroutines ***
60
62
65
66CONTAINS
67
68! **************************************************************************************************
69!> \brief Initialize all the atomic kind radii for a given threshold value.
70!> \param qs_control ...
71!> \param qs_kind_set ...
72!> \date 24.09.2002
73!> \author GT
74!> \version 1.0
75! **************************************************************************************************
76 SUBROUTINE init_interaction_radii(qs_control, qs_kind_set)
77
78 TYPE(qs_control_type), INTENT(IN) :: qs_control
79 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
80
81 CHARACTER(len=*), PARAMETER :: routinen = 'init_interaction_radii'
82
83 INTEGER :: handle, i, iexp_ppl, ikind, ip, iprj_ppnl, j, l, lppl, lppnl, lppsl, lprj, &
84 lprj_ppnl, maxl, n_local, nexp_lpot, nexp_lsd, nexp_ppl, nkind, nloc
85 INTEGER, DIMENSION(0:10) :: npot
86 INTEGER, DIMENSION(1:10) :: nrloc
87 INTEGER, DIMENSION(1:15, 0:10) :: nrpot
88 INTEGER, DIMENSION(:), POINTER :: nct_lpot, nct_lsd, nprj, nprj_ppnl
89 LOGICAL :: ecp_local, ecp_semi_local, llsd, lpot, &
90 paw_atom
91 LOGICAL, DIMENSION(0:5) :: is_nonlocal
92 REAL(kind=dp) :: alpha_core_charge, alpha_ppl, ccore_charge, cerf_ppl, core_charge_radius, &
93 cval, ppl_radius, ppnl_radius, rcprj, zeta
94 REAL(kind=dp), DIMENSION(1:10) :: aloc, bloc
95 REAL(kind=dp), DIMENSION(1:15, 0:10) :: apot, bpot
96 REAL(kind=dp), DIMENSION(:), POINTER :: a_local, a_nonlocal, alpha_lpot, &
97 alpha_lsd, alpha_ppnl, c_local, &
98 cexp_ppl
99 REAL(kind=dp), DIMENSION(:, :), POINTER :: cprj_ppnl, cval_lpot, cval_lsd, rzetprj, &
100 zet
101 REAL(kind=dp), DIMENSION(:, :, :), POINTER :: c_nonlocal
102 TYPE(all_potential_type), POINTER :: all_potential
103 TYPE(cneo_potential_type), POINTER :: cneo_potential
104 TYPE(gth_potential_type), POINTER :: gth_potential
105 TYPE(gto_basis_set_type), POINTER :: aux_basis_set, aux_fit_basis_set, aux_gw_basis, &
106 aux_opt_basis_set, gapw_1c_basis, harris_basis, lri_basis, mao_basis, min_basis_set, &
107 nuc_basis_set, orb_basis_set, p_lri_basis, rhoin_basis, ri_aux_basis_set, ri_basis, &
108 ri_xas_basis, soft_basis, tda_k_basis
109 TYPE(paw_proj_set_type), POINTER :: paw_proj_set
110 TYPE(sgp_potential_type), POINTER :: sgp_potential
111
112 CALL timeset(routinen, handle)
113
114 NULLIFY (all_potential, gth_potential, sgp_potential, cneo_potential)
115 NULLIFY (aux_basis_set, aux_fit_basis_set, aux_gw_basis, tda_k_basis, &
116 harris_basis, lri_basis, mao_basis, orb_basis_set, p_lri_basis, ri_aux_basis_set, &
117 ri_basis, ri_xas_basis, soft_basis, gapw_1c_basis, aux_opt_basis_set, min_basis_set, &
118 nuc_basis_set)
119 NULLIFY (nprj_ppnl, nprj)
120 NULLIFY (alpha_ppnl, cexp_ppl, cprj_ppnl, zet)
121
122 nkind = SIZE(qs_kind_set)
123 nexp_ppl = 0
124
125 DO ikind = 1, nkind
126
127 CALL get_qs_kind(qs_kind_set(ikind), basis_set=orb_basis_set, basis_type="ORB")
128 CALL get_qs_kind(qs_kind_set(ikind), basis_set=aux_basis_set, basis_type="AUX")
129 CALL get_qs_kind(qs_kind_set(ikind), basis_set=min_basis_set, basis_type="MIN")
130 CALL get_qs_kind(qs_kind_set(ikind), basis_set=aux_fit_basis_set, basis_type="AUX_FIT")
131 CALL get_qs_kind(qs_kind_set(ikind), basis_set=aux_opt_basis_set, basis_type="AUX_OPT")
132 CALL get_qs_kind(qs_kind_set(ikind), basis_set=lri_basis, basis_type="LRI_AUX")
133 CALL get_qs_kind(qs_kind_set(ikind), basis_set=p_lri_basis, basis_type="P_LRI_AUX")
134 CALL get_qs_kind(qs_kind_set(ikind), basis_set=ri_basis, basis_type="RI_HXC")
135 CALL get_qs_kind(qs_kind_set(ikind), basis_set=ri_aux_basis_set, basis_type="RI_AUX")
136 CALL get_qs_kind(qs_kind_set(ikind), basis_set=mao_basis, basis_type="MAO")
137 CALL get_qs_kind(qs_kind_set(ikind), basis_set=harris_basis, basis_type="HARRIS")
138 CALL get_qs_kind(qs_kind_set(ikind), basis_set=aux_gw_basis, basis_type="AUX_GW")
139 CALL get_qs_kind(qs_kind_set(ikind), basis_set=ri_xas_basis, basis_type="RI_XAS")
140 CALL get_qs_kind(qs_kind_set(ikind), basis_set=soft_basis, basis_type="ORB_SOFT")
141 CALL get_qs_kind(qs_kind_set(ikind), basis_set=gapw_1c_basis, basis_type="GAPW_1C")
142 CALL get_qs_kind(qs_kind_set(ikind), basis_set=tda_k_basis, basis_type="TDA_HFX")
143 CALL get_qs_kind(qs_kind_set(ikind), basis_set=rhoin_basis, basis_type="RHOIN")
144 CALL get_qs_kind(qs_kind_set(ikind), basis_set=nuc_basis_set, basis_type="NUC")
145 CALL get_qs_kind(qs_kind_set(ikind), &
146 paw_proj_set=paw_proj_set, &
147 paw_atom=paw_atom, &
148 all_potential=all_potential, &
149 gth_potential=gth_potential, &
150 sgp_potential=sgp_potential, &
151 cneo_potential=cneo_potential)
152
153 ! Calculate the orbital basis function radii ***
154 ! For ALL electron this has to come before the calculation of the
155 ! radii used to calculate the 3c lists.
156 ! Indeed, there the radii of the GTO primitives are needed
157 IF (ASSOCIATED(orb_basis_set)) THEN
158 CALL init_interaction_radii_orb_basis(orb_basis_set, qs_control%eps_pgf_orb, &
159 qs_control%eps_kg_orb)
160 END IF
161
162 IF (ASSOCIATED(all_potential)) THEN
163
164 CALL get_potential(potential=all_potential, &
165 alpha_core_charge=alpha_core_charge, &
166 ccore_charge=ccore_charge)
167
168 ! Calculate the radius of the core charge distribution
169
170 core_charge_radius = exp_radius(0, alpha_core_charge, &
171 qs_control%eps_core_charge, &
172 ccore_charge)
173
174 CALL set_potential(potential=all_potential, &
175 core_charge_radius=core_charge_radius)
176
177 ELSE IF (ASSOCIATED(gth_potential)) THEN
178
179 CALL get_potential(potential=gth_potential, &
180 alpha_core_charge=alpha_core_charge, &
181 ccore_charge=ccore_charge, &
182 alpha_ppl=alpha_ppl, &
183 cerf_ppl=cerf_ppl, &
184 nexp_ppl=nexp_ppl, &
185 cexp_ppl=cexp_ppl, &
186 lpot_present=lpot, &
187 lsd_present=llsd, &
188 lppnl=lppnl, &
189 alpha_ppnl=alpha_ppnl, &
190 nprj_ppnl=nprj_ppnl, &
191 cprj_ppnl=cprj_ppnl)
192
193 ! Calculate the radius of the core charge distribution
194 core_charge_radius = exp_radius(0, alpha_core_charge, &
195 qs_control%eps_core_charge, &
196 ccore_charge)
197
198 ! Calculate the radii of the local part
199 ! of the Goedecker pseudopotential (GTH)
200 ppl_radius = exp_radius(0, alpha_ppl, qs_control%eps_ppl, cerf_ppl)
201
202 DO iexp_ppl = 1, nexp_ppl
203 lppl = 2*(iexp_ppl - 1)
204 ppl_radius = max(ppl_radius, &
205 exp_radius(lppl, alpha_ppl, &
206 qs_control%eps_ppl, &
207 cexp_ppl(iexp_ppl), &
208 rlow=ppl_radius))
209 END DO
210 IF (lpot) THEN
211 ! extended local potential
212 CALL get_potential(potential=gth_potential, &
213 nexp_lpot=nexp_lpot, &
214 alpha_lpot=alpha_lpot, &
215 nct_lpot=nct_lpot, &
216 cval_lpot=cval_lpot)
217 DO j = 1, nexp_lpot
218 DO i = 1, nct_lpot(j)
219 lppl = 2*(i - 1)
220 ppl_radius = max(ppl_radius, &
221 exp_radius(lppl, alpha_lpot(j), qs_control%eps_ppl, &
222 cval_lpot(i, j), rlow=ppl_radius))
223 END DO
224 END DO
225 END IF
226 IF (llsd) THEN
227 ! lsd local potential
228 CALL get_potential(potential=gth_potential, &
229 nexp_lsd=nexp_lsd, &
230 alpha_lsd=alpha_lsd, &
231 nct_lsd=nct_lsd, &
232 cval_lsd=cval_lsd)
233 DO j = 1, nexp_lsd
234 DO i = 1, nct_lsd(j)
235 lppl = 2*(i - 1)
236 ppl_radius = max(ppl_radius, &
237 exp_radius(lppl, alpha_lsd(j), qs_control%eps_ppl, &
238 cval_lsd(i, j), rlow=ppl_radius))
239 END DO
240 END DO
241 END IF
242
243 ! Calculate the radii of the non-local part
244 ! of the Goedecker pseudopotential (GTH)
245 ppnl_radius = 0.0_dp
246 DO l = 0, lppnl
247 DO iprj_ppnl = 1, nprj_ppnl(l)
248 lprj_ppnl = l + 2*(iprj_ppnl - 1)
249 ppnl_radius = max(ppnl_radius, &
250 exp_radius(lprj_ppnl, alpha_ppnl(l), &
251 qs_control%eps_pgf_orb, &
252 cprj_ppnl(iprj_ppnl, l), &
253 rlow=ppnl_radius))
254 END DO
255 END DO
256 CALL set_potential(potential=gth_potential, &
257 core_charge_radius=core_charge_radius, &
258 ppl_radius=ppl_radius, &
259 ppnl_radius=ppnl_radius)
260
261 ELSE IF (ASSOCIATED(sgp_potential)) THEN
262
263 ! Calculate the radius of the core charge distribution
264 CALL get_potential(potential=sgp_potential, &
265 alpha_core_charge=alpha_core_charge, &
266 ccore_charge=ccore_charge)
267 core_charge_radius = exp_radius(0, alpha_core_charge, &
268 qs_control%eps_core_charge, &
269 ccore_charge)
270 ! Calculate the radii of the local part
271 ppl_radius = core_charge_radius
272 CALL get_potential(potential=sgp_potential, ecp_local=ecp_local)
273 IF (ALLOCATED(sgp_potential%upf_data)) THEN
274 IF (sgp_potential%upf_local%radial%degree < 0 .OR. &
275 sgp_potential%upf_local%eps /= qs_control%eps_ppl) THEN
276 CALL upf_local_init(sgp_potential%upf_data, alpha_core_charge, &
277 max(0, sgp_potential%upf_local%radial%degree), &
278 sgp_potential%upf_local, qs_control%eps_ppl)
279 END IF
280 ppl_radius = max(ppl_radius, sgp_potential%upf_local%radial%radius)
281 ELSE IF (ecp_local) THEN
282 CALL get_potential(potential=sgp_potential, nloc=nloc, nrloc=nrloc, aloc=aloc, bloc=bloc)
283 DO i = 1, nloc
284 lppl = max(0, nrloc(i) - 2)
285 ppl_radius = max(ppl_radius, &
286 exp_radius(lppl, bloc(i), qs_control%eps_ppl, aloc(i), rlow=ppl_radius))
287 END DO
288 ELSE
289 CALL get_potential(potential=sgp_potential, n_local=n_local, a_local=a_local, c_local=c_local)
290 DO i = 1, n_local
291 ppl_radius = max(ppl_radius, &
292 exp_radius(0, a_local(i), qs_control%eps_ppl, c_local(i), rlow=ppl_radius))
293 END DO
294 END IF
295 ! Calculate the radii of the semi-local part
296 CALL get_potential(potential=sgp_potential, ecp_semi_local=ecp_semi_local)
297 IF (ecp_semi_local) THEN
298 CALL get_potential(potential=sgp_potential, sl_lmax=lppsl, npot=npot, nrpot=nrpot, &
299 apot=apot, bpot=bpot)
300 DO l = 0, lppsl
301 DO i = 1, npot(l)
302 lppl = max(0, nrpot(i, l) - 2)
303 ppl_radius = max(ppl_radius, &
304 exp_radius(lppl, bpot(i, l), qs_control%eps_ppl, apot(i, l), &
305 rlow=ppl_radius))
306 END DO
307 END DO
308 END IF
309 ! Calculate the radii of the non-local part
310 ppnl_radius = 0.0_dp
311 CALL get_potential(potential=sgp_potential, lmax=lppnl, n_nonlocal=nloc)
312 CALL get_potential(potential=sgp_potential, is_nonlocal=is_nonlocal, &
313 a_nonlocal=a_nonlocal, c_nonlocal=c_nonlocal)
314 IF (ALLOCATED(sgp_potential%upf_data)) THEN
315 ppnl_radius = sgp_potential%upf_projectors%radius
316 ELSE
317 DO l = 0, lppnl
318 IF (is_nonlocal(l)) THEN
319 DO i = 1, nloc
320 ! SGP coefficients contain 1/sqrt(4*pi). The addition
321 ! theorem bounds every angular component by sqrt(2*l+1).
322 cval = sqrt(real(2*l + 1, dp))*maxval(abs(c_nonlocal(i, :, l)))
323 ppnl_radius = max(ppnl_radius, &
324 exp_radius(l, a_nonlocal(i), qs_control%eps_pgf_orb, cval, rlow=ppnl_radius))
325 END DO
326 END IF
327 END DO
328 END IF
329 CALL set_potential(potential=sgp_potential, &
330 core_charge_radius=core_charge_radius, &
331 ppl_radius=ppl_radius, &
332 ppnl_radius=ppnl_radius)
333
334 ELSE IF (ASSOCIATED(cneo_potential)) THEN
335
336 IF (ASSOCIATED(nuc_basis_set)) THEN
337 CALL init_interaction_radii_orb_basis(nuc_basis_set, qs_control%eps_pgf_orb/ &
338 sqrt(cneo_potential%zeff))
339 END IF
340
341 END IF
342
343 ! Calculate the aux fit orbital basis function radii
344 IF (ASSOCIATED(aux_fit_basis_set)) THEN
345 CALL init_interaction_radii_orb_basis(aux_fit_basis_set, qs_control%eps_pgf_orb, &
346 qs_control%eps_kg_orb)
347 END IF
348
349 ! Calculate the aux opt orbital basis function radii
350 IF (ASSOCIATED(aux_opt_basis_set)) THEN
351 CALL init_interaction_radii_orb_basis(aux_opt_basis_set, qs_control%eps_pgf_orb)
352 END IF
353
354 ! Calculate the minimal basis function radii
355 IF (ASSOCIATED(min_basis_set)) THEN
356 CALL init_interaction_radii_orb_basis(min_basis_set, qs_control%eps_pgf_orb)
357 END IF
358
359 ! Calculate the aux orbital basis function radii
360 IF (ASSOCIATED(aux_basis_set)) THEN
361 CALL init_interaction_radii_orb_basis(aux_basis_set, qs_control%eps_pgf_orb)
362 END IF
363
364 ! Calculate the RI aux orbital basis function radii
365 IF (ASSOCIATED(ri_aux_basis_set)) THEN
366 CALL init_interaction_radii_orb_basis(ri_aux_basis_set, qs_control%eps_pgf_orb)
367 END IF
368
369 ! Calculate the MAO orbital basis function radii
370 IF (ASSOCIATED(mao_basis)) THEN
371 CALL init_interaction_radii_orb_basis(mao_basis, qs_control%eps_pgf_orb)
372 END IF
373
374 ! Calculate the HARRIS orbital basis function radii
375 IF (ASSOCIATED(harris_basis)) THEN
376 CALL init_interaction_radii_orb_basis(harris_basis, qs_control%eps_pgf_orb)
377 END IF
378
379 ! Calculate the AUX_GW orbital basis function radii
380 IF (ASSOCIATED(aux_gw_basis)) THEN
381 CALL init_interaction_radii_orb_basis(aux_gw_basis, qs_control%eps_pgf_orb)
382 END IF
383
384 ! Calculate the soft orbital basis function radii
385 IF (ASSOCIATED(soft_basis)) THEN
386 CALL init_interaction_radii_orb_basis(soft_basis, qs_control%eps_pgf_orb)
387 END IF
388
389 ! Calculate the GAPW 1 center basis function radii
390 IF (ASSOCIATED(gapw_1c_basis)) THEN
391 CALL init_interaction_radii_orb_basis(gapw_1c_basis, qs_control%eps_pgf_orb)
392 END IF
393
394 ! Calculate the HFX SR basis for TDA
395 IF (ASSOCIATED(tda_k_basis)) THEN
396 CALL init_interaction_radii_orb_basis(tda_k_basis, qs_control%eps_pgf_orb)
397 END IF
398
399 ! Calculate the lri basis function radii
400 IF (ASSOCIATED(lri_basis)) THEN
401 CALL init_interaction_radii_orb_basis(lri_basis, qs_control%eps_pgf_orb)
402 END IF
403 IF (ASSOCIATED(ri_basis)) THEN
404 CALL init_interaction_radii_orb_basis(ri_basis, qs_control%eps_pgf_orb)
405 END IF
406 IF (ASSOCIATED(ri_xas_basis)) THEN
407 CALL init_interaction_radii_orb_basis(ri_xas_basis, qs_control%eps_pgf_orb)
408 END IF
409 IF (ASSOCIATED(p_lri_basis)) THEN
410 CALL init_interaction_radii_orb_basis(p_lri_basis, qs_control%eps_pgf_orb)
411 END IF
412
413 ! Calculate the density input basis function radii
414 IF (ASSOCIATED(rhoin_basis)) THEN
415 CALL init_interaction_radii_orb_basis(rhoin_basis, qs_control%eps_pgf_orb)
416 END IF
417
418 ! Calculate the paw projector basis function radii
419 IF (ASSOCIATED(paw_proj_set)) THEN
420 CALL get_paw_proj_set(paw_proj_set=paw_proj_set, &
421 nprj=nprj, &
422 maxl=maxl, &
423 zetprj=zet, &
424 rzetprj=rzetprj)
425 rcprj = 0.0_dp
426 rzetprj = 0.0_dp
427 DO lprj = 0, maxl
428 DO ip = 1, nprj(lprj)
429 zeta = zet(ip, lprj)
430 rzetprj(ip, lprj) = max(rzetprj(ip, lprj), &
431 exp_radius(lprj, zeta, qs_control%eps_pgf_orb, &
432 1.0_dp, rlow=rzetprj(ip, lprj)))
433 rcprj = max(rcprj, rzetprj(ip, lprj))
434 END DO
435 END DO
436 CALL set_paw_proj_set(paw_proj_set=paw_proj_set, &
437 rzetprj=rzetprj, &
438 rcprj=rcprj)
439
440 END IF
441 END DO ! ikind
442
443 CALL timestop(handle)
444
445 END SUBROUTINE init_interaction_radii
446
447! **************************************************************************************************
448!> \brief ...
449!> \param orb_basis_set ...
450!> \param eps_pgf_orb ...
451!> \param eps_pgf_short ...
452! **************************************************************************************************
453 SUBROUTINE init_interaction_radii_orb_basis(orb_basis_set, eps_pgf_orb, eps_pgf_short)
454
455 TYPE(gto_basis_set_type), POINTER :: orb_basis_set
456 REAL(kind=dp), INTENT(IN) :: eps_pgf_orb
457 REAL(kind=dp), INTENT(IN), OPTIONAL :: eps_pgf_short
458
459 INTEGER :: ipgf, iset, ishell, l, nset
460 INTEGER, DIMENSION(:), POINTER :: npgf, nshell
461 INTEGER, DIMENSION(:, :), POINTER :: lshell
462 REAL(kind=dp) :: eps_short, gcca, kind_radius, &
463 short_radius, zeta
464 REAL(kind=dp), DIMENSION(:), POINTER :: set_radius
465 REAL(kind=dp), DIMENSION(:, :), POINTER :: pgf_radius, zet
466 REAL(kind=dp), DIMENSION(:, :, :), POINTER :: gcc
467
468 IF (ASSOCIATED(orb_basis_set)) THEN
469
470 IF (PRESENT(eps_pgf_short)) THEN
471 eps_short = eps_pgf_short
472 ELSE
473 eps_short = eps_pgf_orb
474 END IF
475
476 CALL get_gto_basis_set(gto_basis_set=orb_basis_set, &
477 nset=nset, &
478 nshell=nshell, &
479 npgf=npgf, &
480 l=lshell, &
481 zet=zet, &
482 gcc=gcc, &
483 pgf_radius=pgf_radius, &
484 set_radius=set_radius)
485
486 kind_radius = 0.0_dp
487 short_radius = 0.0_dp
488
489 DO iset = 1, nset
490 set_radius(iset) = 0.0_dp
491 DO ipgf = 1, npgf(iset)
492 pgf_radius(ipgf, iset) = 0.0_dp
493 DO ishell = 1, nshell(iset)
494 l = lshell(ishell, iset)
495 gcca = gcc(ipgf, ishell, iset)
496 zeta = zet(ipgf, iset)
497 pgf_radius(ipgf, iset) = max(pgf_radius(ipgf, iset), &
498 exp_radius(l, zeta, eps_pgf_orb, gcca, &
499 rlow=pgf_radius(ipgf, iset)))
500 short_radius = max(short_radius, &
501 exp_radius(l, zeta, eps_short, gcca, &
502 rlow=short_radius))
503 END DO
504 set_radius(iset) = max(set_radius(iset), pgf_radius(ipgf, iset))
505 END DO
506 kind_radius = max(kind_radius, set_radius(iset))
507 END DO
508
509 CALL set_gto_basis_set(gto_basis_set=orb_basis_set, &
510 pgf_radius=pgf_radius, &
511 set_radius=set_radius, &
512 kind_radius=kind_radius, &
513 short_kind_radius=short_radius)
514
515 END IF
516
518
519! **************************************************************************************************
520!> \brief ...
521!> \param se_control ...
522!> \param atomic_kind_set ...
523!> \param qs_kind_set ...
524!> \param subsys_section ...
525! **************************************************************************************************
526 SUBROUTINE init_se_nlradius(se_control, atomic_kind_set, qs_kind_set, &
527 subsys_section)
528
529 TYPE(semi_empirical_control_type), POINTER :: se_control
530 TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
531 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
532 TYPE(section_vals_type), POINTER :: subsys_section
533
534 INTEGER :: ikind, nkind
535 REAL(kind=dp) :: kind_radius
536 TYPE(gto_basis_set_type), POINTER :: orb_basis_set
537 TYPE(qs_kind_type), POINTER :: qs_kind
538
539 NULLIFY (orb_basis_set, qs_kind)
540
541 nkind = SIZE(qs_kind_set)
542
543 DO ikind = 1, nkind
544
545 qs_kind => qs_kind_set(ikind)
546
547 CALL get_qs_kind(qs_kind=qs_kind, basis_set=orb_basis_set)
548
549 IF (ASSOCIATED(orb_basis_set)) THEN
550
551 CALL get_gto_basis_set(gto_basis_set=orb_basis_set, &
552 kind_radius=kind_radius)
553
554 kind_radius = max(kind_radius, se_control%cutoff_exc)
555
556 CALL set_gto_basis_set(gto_basis_set=orb_basis_set, &
557 kind_radius=kind_radius)
558
559 END IF
560
561 END DO
562
563 CALL write_kind_radii(atomic_kind_set, qs_kind_set, subsys_section)
564
565 END SUBROUTINE init_se_nlradius
566
567! **************************************************************************************************
568!> \brief ...
569!> \param atomic_kind_set ...
570!> \param qs_kind_set ...
571!> \param subsys_section ...
572! **************************************************************************************************
573 SUBROUTINE write_kind_radii(atomic_kind_set, qs_kind_set, subsys_section)
574
575 TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
576 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
577 TYPE(section_vals_type), POINTER :: subsys_section
578
579 CHARACTER(LEN=10) :: bas
580 CHARACTER(LEN=default_string_length) :: name, unit_str
581 INTEGER :: ikind, nkind, output_unit
582 REAL(kind=dp) :: conv, kind_radius
583 TYPE(cp_logger_type), POINTER :: logger
584 TYPE(gto_basis_set_type), POINTER :: orb_basis_set
585
586 NULLIFY (logger)
587 logger => cp_get_default_logger()
588 NULLIFY (orb_basis_set)
589 bas = "ORBITAL "
590
591 nkind = SIZE(atomic_kind_set)
592! *** Print the kind radii ***
593 output_unit = cp_print_key_unit_nr(logger, subsys_section, &
594 "PRINT%RADII/KIND_RADII", extension=".Log")
595 CALL section_vals_val_get(subsys_section, "PRINT%RADII%UNIT", c_val=unit_str)
596 conv = cp_unit_from_cp2k(1.0_dp, trim(unit_str))
597 IF (output_unit > 0) THEN
598 WRITE (unit=output_unit, fmt="(/,T2,A,T56,A,T63,A,T75,A)") &
599 "RADII: "//trim(bas)//" BASIS in "//trim(unit_str), &
600 "Kind", "Label", "Radius"
601 DO ikind = 1, nkind
602 CALL get_atomic_kind(atomic_kind_set(ikind), name=name)
603 CALL get_qs_kind(qs_kind_set(ikind), basis_set=orb_basis_set)
604 IF (ASSOCIATED(orb_basis_set)) THEN
605 CALL get_gto_basis_set(gto_basis_set=orb_basis_set, &
606 kind_radius=kind_radius)
607 WRITE (unit=output_unit, fmt="(T55,I5,T63,A5,T69,F12.6)") &
608 ikind, name, kind_radius*conv
609 ELSE
610 WRITE (unit=output_unit, fmt="(T55,I5,T63,A5,T72,A9)") &
611 ikind, name, "no basis"
612 END IF
613 END DO
614 END IF
615 CALL cp_print_key_finished_output(output_unit, logger, subsys_section, &
616 "PRINT%RADII/KIND_RADII")
617
618 END SUBROUTINE write_kind_radii
619
620! **************************************************************************************************
621!> \brief Write the radii of the core charge distributions to the output
622!> unit.
623!> \param atomic_kind_set ...
624!> \param qs_kind_set ...
625!> \param subsys_section ...
626!> \date 15.09.2000
627!> \author MK
628!> \version 1.0
629! **************************************************************************************************
630 SUBROUTINE write_core_charge_radii(atomic_kind_set, qs_kind_set, subsys_section)
631 TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
632 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
633 TYPE(section_vals_type), POINTER :: subsys_section
634
635 CHARACTER(LEN=default_string_length) :: name, unit_str
636 INTEGER :: ikind, nkind, output_unit
637 REAL(kind=dp) :: conv, core_charge_radius
638 TYPE(all_potential_type), POINTER :: all_potential
639 TYPE(cp_logger_type), POINTER :: logger
640 TYPE(gth_potential_type), POINTER :: gth_potential
641
642 NULLIFY (logger)
643 logger => cp_get_default_logger()
644 output_unit = cp_print_key_unit_nr(logger, subsys_section, &
645 "PRINT%RADII/CORE_CHARGE_RADII", extension=".Log")
646 CALL section_vals_val_get(subsys_section, "PRINT%RADII%UNIT", c_val=unit_str)
647 conv = cp_unit_from_cp2k(1.0_dp, trim(unit_str))
648 IF (output_unit > 0) THEN
649 nkind = SIZE(atomic_kind_set)
650 WRITE (unit=output_unit, fmt="(/,T2,A,T56,A,T63,A,T75,A)") &
651 "RADII: CORE CHARGE DISTRIBUTIONS in "// &
652 trim(unit_str), "Kind", "Label", "Radius"
653 DO ikind = 1, nkind
654 CALL get_atomic_kind(atomic_kind_set(ikind), name=name)
655 CALL get_qs_kind(qs_kind_set(ikind), &
656 all_potential=all_potential, gth_potential=gth_potential)
657
658 IF (ASSOCIATED(all_potential)) THEN
659 CALL get_potential(potential=all_potential, &
660 core_charge_radius=core_charge_radius)
661 WRITE (unit=output_unit, fmt="(T55,I5,T63,A5,T69,F12.6)") &
662 ikind, name, core_charge_radius*conv
663 ELSE IF (ASSOCIATED(gth_potential)) THEN
664 CALL get_potential(potential=gth_potential, &
665 core_charge_radius=core_charge_radius)
666 WRITE (unit=output_unit, fmt="(T55,I5,T63,A5,T69,F12.6)") &
667 ikind, name, core_charge_radius*conv
668 END IF
669 END DO
670 END IF
671 CALL cp_print_key_finished_output(output_unit, logger, subsys_section, &
672 "PRINT%RADII/CORE_CHARGE_RADII")
673 END SUBROUTINE write_core_charge_radii
674
675! **************************************************************************************************
676!> \brief Write the orbital basis function radii to the output unit.
677!> \param basis ...
678!> \param atomic_kind_set ...
679!> \param qs_kind_set ...
680!> \param subsys_section ...
681!> \date 05.06.2000
682!> \author MK
683!> \version 1.0
684! **************************************************************************************************
685 SUBROUTINE write_pgf_orb_radii(basis, atomic_kind_set, qs_kind_set, subsys_section)
686
687 CHARACTER(len=*) :: basis
688 TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
689 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
690 TYPE(section_vals_type), POINTER :: subsys_section
691
692 CHARACTER(LEN=10) :: bas
693 CHARACTER(LEN=default_string_length) :: name, unit_str
694 INTEGER :: ikind, ipgf, iset, nkind, nset, &
695 output_unit
696 INTEGER, DIMENSION(:), POINTER :: npgf
697 REAL(kind=dp) :: conv, kind_radius, short_kind_radius
698 REAL(kind=dp), DIMENSION(:), POINTER :: set_radius
699 REAL(kind=dp), DIMENSION(:, :), POINTER :: pgf_radius
700 TYPE(cp_logger_type), POINTER :: logger
701 TYPE(gto_basis_set_type), POINTER :: aux_basis_set, lri_basis_set, &
702 nuc_basis_set, orb_basis_set
703
704 NULLIFY (logger)
705 logger => cp_get_default_logger()
706 NULLIFY (aux_basis_set, orb_basis_set, lri_basis_set, nuc_basis_set)
707 bas = " "
708 bas(1:3) = basis(1:3)
709 CALL uppercase(bas)
710 IF (bas(1:3) == "ORB") THEN
711 bas = "ORBITAL "
712 ELSE IF (bas(1:3) == "AUX") THEN
713 bas = "AUXILLIARY"
714 ELSE IF (bas(1:3) == "LRI") THEN
715 bas = "LOCAL RI"
716 ELSE IF (bas(1:3) == "NUC") THEN
717 bas = "NUCLEAR "
718 ELSE
719 cpabort("Undefined basis set type")
720 END IF
721
722 nkind = SIZE(qs_kind_set)
723
724! *** Print the kind radii ***
725 output_unit = cp_print_key_unit_nr(logger, subsys_section, &
726 "PRINT%RADII/KIND_RADII", extension=".Log")
727 CALL section_vals_val_get(subsys_section, "PRINT%RADII%UNIT", c_val=unit_str)
728 conv = cp_unit_from_cp2k(1.0_dp, trim(unit_str))
729 IF (output_unit > 0) THEN
730 WRITE (unit=output_unit, fmt="(/,T2,A,T46,A,T53,A,T63,A,T71,A)") &
731 "RADII: "//trim(bas)//" BASIS in "//trim(unit_str), &
732 "Kind", "Label", "Radius", "OCE Radius"
733 DO ikind = 1, nkind
734 CALL get_atomic_kind(atomic_kind_set(ikind), name=name)
735 IF (bas(1:3) == "ORB") THEN
736 CALL get_qs_kind(qs_kind_set(ikind), basis_set=orb_basis_set)
737 IF (ASSOCIATED(orb_basis_set)) THEN
738 CALL get_gto_basis_set(gto_basis_set=orb_basis_set, &
739 kind_radius=kind_radius, &
740 short_kind_radius=short_kind_radius)
741 END IF
742 ELSE IF (bas(1:3) == "AUX") THEN
743 CALL get_qs_kind(qs_kind_set(ikind), basis_set=aux_basis_set, basis_type="AUX")
744 IF (ASSOCIATED(aux_basis_set)) THEN
745 CALL get_gto_basis_set(gto_basis_set=aux_basis_set, &
746 kind_radius=kind_radius)
747 END IF
748 ELSE IF (bas(1:3) == "LOC") THEN
749 CALL get_qs_kind(qs_kind_set(ikind), basis_set=lri_basis_set, basis_type="LRI_AUX")
750 IF (ASSOCIATED(lri_basis_set)) THEN
751 CALL get_gto_basis_set(gto_basis_set=lri_basis_set, &
752 kind_radius=kind_radius)
753 END IF
754 ELSE IF (bas(1:3) == "NUC") THEN
755 CALL get_qs_kind(qs_kind_set(ikind), basis_set=nuc_basis_set, basis_type="NUC")
756 IF (ASSOCIATED(nuc_basis_set)) THEN
757 CALL get_gto_basis_set(gto_basis_set=nuc_basis_set, &
758 kind_radius=kind_radius, &
759 short_kind_radius=short_kind_radius)
760 END IF
761 ELSE
762 cpabort("Undefined basis set type")
763 END IF
764 IF (ASSOCIATED(aux_basis_set) .OR. ASSOCIATED(orb_basis_set) .OR. &
765 ASSOCIATED(nuc_basis_set)) THEN
766 WRITE (unit=output_unit, fmt="(T45,I5,T53,A5,T57,F12.6,T69,F12.6)") &
767 ikind, name, kind_radius*conv, short_kind_radius*conv
768 ELSE
769 WRITE (unit=output_unit, fmt="(T55,I5,T63,A5,T72,A9)") &
770 ikind, name, "no basis"
771 END IF
772 END DO
773 END IF
774 CALL cp_print_key_finished_output(output_unit, logger, subsys_section, &
775 "PRINT%RADII/KIND_RADII")
776
777! *** Print the shell set radii ***
778 output_unit = cp_print_key_unit_nr(logger, subsys_section, &
779 "PRINT%RADII/SET_RADII", extension=".Log")
780 IF (output_unit > 0) THEN
781 WRITE (unit=output_unit, fmt="(/,T2,A,T51,A,T57,A,T65,A,T75,A)") &
782 "RADII: SHELL SETS OF "//trim(bas)//" BASIS in "// &
783 trim(unit_str), "Kind", "Label", "Set", "Radius"
784 DO ikind = 1, nkind
785 CALL get_atomic_kind(atomic_kind_set(ikind), name=name)
786 IF (bas(1:3) == "ORB") THEN
787 CALL get_qs_kind(qs_kind_set(ikind), basis_set=orb_basis_set)
788 IF (ASSOCIATED(orb_basis_set)) THEN
789 CALL get_gto_basis_set(gto_basis_set=orb_basis_set, &
790 nset=nset, &
791 set_radius=set_radius)
792 END IF
793 ELSE IF (bas(1:3) == "AUX") THEN
794 CALL get_qs_kind(qs_kind_set(ikind), basis_set=aux_basis_set, basis_type="AUX")
795 IF (ASSOCIATED(aux_basis_set)) THEN
796 CALL get_gto_basis_set(gto_basis_set=aux_basis_set, &
797 nset=nset, &
798 set_radius=set_radius)
799 END IF
800 ELSE IF (bas(1:3) == "LOC") THEN
801 CALL get_qs_kind(qs_kind_set(ikind), basis_set=lri_basis_set, basis_type="LRI_AUX")
802 IF (ASSOCIATED(lri_basis_set)) THEN
803 CALL get_gto_basis_set(gto_basis_set=lri_basis_set, &
804 nset=nset, &
805 set_radius=set_radius)
806 END IF
807 ELSE IF (bas(1:3) == "NUC") THEN
808 CALL get_qs_kind(qs_kind_set(ikind), basis_set=nuc_basis_set, basis_type="NUC")
809 IF (ASSOCIATED(nuc_basis_set)) THEN
810 CALL get_gto_basis_set(gto_basis_set=nuc_basis_set, &
811 nset=nset, &
812 set_radius=set_radius)
813 END IF
814 ELSE
815 cpabort("Undefined basis set type")
816 END IF
817 IF (ASSOCIATED(aux_basis_set) .OR. ASSOCIATED(orb_basis_set) .OR. &
818 ASSOCIATED(nuc_basis_set)) THEN
819 WRITE (unit=output_unit, fmt="(T50,I5,T57,A5,(T63,I5,T69,F12.6))") &
820 ikind, name, (iset, set_radius(iset)*conv, iset=1, nset)
821 ELSE
822 WRITE (unit=output_unit, fmt="(T50,I5,T58,A5,T73,A8)") &
823 ikind, name, "no basis"
824 END IF
825 END DO
826 END IF
827 CALL cp_print_key_finished_output(output_unit, logger, subsys_section, &
828 "PRINT%RADII/SET_RADII")
829! *** Print the primitive Gaussian function radii ***
830 output_unit = cp_print_key_unit_nr(logger, subsys_section, &
831 "PRINT%RADII/PGF_RADII", extension=".Log")
832 IF (output_unit > 0) THEN
833 WRITE (unit=output_unit, fmt="(/,T2,A,T51,A,T57,A,T65,A,T75,A)") &
834 "RADII: PRIMITIVE GAUSSIANS OF "//trim(bas)//" BASIS in "// &
835 trim(unit_str), "Kind", "Label", "Set", "Radius"
836 DO ikind = 1, nkind
837 CALL get_atomic_kind(atomic_kind_set(ikind), name=name)
838 IF (bas(1:3) == "ORB") THEN
839 CALL get_qs_kind(qs_kind_set(ikind), basis_set=orb_basis_set)
840 IF (ASSOCIATED(orb_basis_set)) THEN
841 CALL get_gto_basis_set(gto_basis_set=orb_basis_set, &
842 nset=nset, &
843 npgf=npgf, &
844 pgf_radius=pgf_radius)
845 END IF
846 ELSE IF (bas(1:3) == "AUX") THEN
847 CALL get_qs_kind(qs_kind_set(ikind), basis_set=aux_basis_set, basis_type="AUX")
848 IF (ASSOCIATED(aux_basis_set)) THEN
849 CALL get_gto_basis_set(gto_basis_set=aux_basis_set, &
850 nset=nset, &
851 npgf=npgf, &
852 pgf_radius=pgf_radius)
853 END IF
854 ELSE IF (bas(1:3) == "LOC") THEN
855 CALL get_qs_kind(qs_kind_set(ikind), basis_set=lri_basis_set, basis_type="LRI_AUX")
856 IF (ASSOCIATED(lri_basis_set)) THEN
857 CALL get_gto_basis_set(gto_basis_set=lri_basis_set, &
858 nset=nset, &
859 npgf=npgf, &
860 pgf_radius=pgf_radius)
861 END IF
862 ELSE IF (bas(1:3) == "NUC") THEN
863 CALL get_qs_kind(qs_kind_set(ikind), basis_set=nuc_basis_set, basis_type="NUC")
864 IF (ASSOCIATED(nuc_basis_set)) THEN
865 CALL get_gto_basis_set(gto_basis_set=nuc_basis_set, &
866 nset=nset, &
867 npgf=npgf, &
868 pgf_radius=pgf_radius)
869 END IF
870 ELSE
871 cpabort("Undefined basis type")
872 END IF
873 IF (ASSOCIATED(aux_basis_set) .OR. ASSOCIATED(orb_basis_set) .OR. &
874 ASSOCIATED(lri_basis_set) .OR. ASSOCIATED(nuc_basis_set)) THEN
875 DO iset = 1, nset
876 WRITE (unit=output_unit, fmt="(T50,I5,T57,A5,T63,I5,(T69,F12.6))") &
877 ikind, name, iset, &
878 (pgf_radius(ipgf, iset)*conv, ipgf=1, npgf(iset))
879 END DO
880 ELSE
881 WRITE (unit=output_unit, fmt="(T50,I5,T58,A5,T73,A8)") &
882 ikind, name, "no basis"
883 END IF
884 END DO
885 END IF
886 CALL cp_print_key_finished_output(output_unit, logger, subsys_section, &
887 "PRINT%RADII/PGF_RADII")
888 END SUBROUTINE write_pgf_orb_radii
889
890! **************************************************************************************************
891!> \brief Write the radii of the exponential functions of the Goedecker
892!> pseudopotential (GTH, local part) to the logical unit number
893!> "output_unit".
894!> \param atomic_kind_set ...
895!> \param qs_kind_set ...
896!> \param subsys_section ...
897!> \date 06.10.2000
898!> \author MK
899!> \version 1.0
900! **************************************************************************************************
901 SUBROUTINE write_ppl_radii(atomic_kind_set, qs_kind_set, subsys_section)
902
903 TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
904 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
905 TYPE(section_vals_type), POINTER :: subsys_section
906
907 CHARACTER(LEN=default_string_length) :: name, unit_str
908 INTEGER :: ikind, nkind, output_unit
909 REAL(kind=dp) :: conv, ppl_radius
910 TYPE(cp_logger_type), POINTER :: logger
911 TYPE(gth_potential_type), POINTER :: gth_potential
912
913 NULLIFY (logger)
914 logger => cp_get_default_logger()
915 output_unit = cp_print_key_unit_nr(logger, subsys_section, &
916 "PRINT%RADII/PPL_RADII", extension=".Log")
917 CALL section_vals_val_get(subsys_section, "PRINT%RADII%UNIT", c_val=unit_str)
918 conv = cp_unit_from_cp2k(1.0_dp, trim(unit_str))
919 IF (output_unit > 0) THEN
920 nkind = SIZE(atomic_kind_set)
921 WRITE (unit=output_unit, fmt="(/,T2,A,T56,A,T63,A,T75,A)") &
922 "RADII: LOCAL PART OF GTH/ELP PP in "// &
923 trim(unit_str), "Kind", "Label", "Radius"
924 DO ikind = 1, nkind
925 CALL get_atomic_kind(atomic_kind_set(ikind), name=name)
926 CALL get_qs_kind(qs_kind_set(ikind), gth_potential=gth_potential)
927 IF (ASSOCIATED(gth_potential)) THEN
928 CALL get_potential(potential=gth_potential, &
929 ppl_radius=ppl_radius)
930 WRITE (unit=output_unit, fmt="(T55,I5,T63,A5,T69,F12.6)") &
931 ikind, name, ppl_radius*conv
932 END IF
933 END DO
934 END IF
935 CALL cp_print_key_finished_output(output_unit, logger, subsys_section, &
936 "PRINT%RADII/PPL_RADII")
937 END SUBROUTINE write_ppl_radii
938
939! **************************************************************************************************
940!> \brief Write the radii of the projector functions of the Goedecker
941!> pseudopotential (GTH, non-local part) to the logical unit number
942!> "output_unit".
943!> \param atomic_kind_set ...
944!> \param qs_kind_set ...
945!> \param subsys_section ...
946!> \date 06.10.2000
947!> \author MK
948!> \version 1.0
949! **************************************************************************************************
950 SUBROUTINE write_ppnl_radii(atomic_kind_set, qs_kind_set, subsys_section)
951
952 TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
953 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
954 TYPE(section_vals_type), POINTER :: subsys_section
955
956 CHARACTER(LEN=default_string_length) :: name, unit_str
957 INTEGER :: ikind, nkind, output_unit
958 REAL(kind=dp) :: conv, ppnl_radius
959 TYPE(cp_logger_type), POINTER :: logger
960 TYPE(gth_potential_type), POINTER :: gth_potential
961
962 NULLIFY (logger)
963 logger => cp_get_default_logger()
964 output_unit = cp_print_key_unit_nr(logger, subsys_section, &
965 "PRINT%RADII/PPNL_RADII", extension=".Log")
966 CALL section_vals_val_get(subsys_section, "PRINT%RADII%UNIT", c_val=unit_str)
967 conv = cp_unit_from_cp2k(1.0_dp, trim(unit_str))
968 IF (output_unit > 0) THEN
969 nkind = SIZE(atomic_kind_set)
970 WRITE (unit=output_unit, fmt="(/,T2,A,T56,A,T63,A,T75,A)") &
971 "RADII: NON-LOCAL PART OF GTH PP in "// &
972 trim(unit_str), "Kind", "Label", "Radius"
973 DO ikind = 1, nkind
974 CALL get_atomic_kind(atomic_kind_set(ikind), name=name)
975 CALL get_qs_kind(qs_kind_set(ikind), gth_potential=gth_potential)
976 IF (ASSOCIATED(gth_potential)) THEN
977 CALL get_potential(potential=gth_potential, &
978 ppnl_radius=ppnl_radius)
979 WRITE (unit=output_unit, fmt="(T55,I5,T63,A5,T69,F12.6)") &
980 ikind, name, ppnl_radius*conv
981 END IF
982 END DO
983 END IF
984 CALL cp_print_key_finished_output(output_unit, logger, subsys_section, &
985 "PRINT%RADII/PPNL_RADII")
986 END SUBROUTINE write_ppnl_radii
987
988! **************************************************************************************************
989!> \brief Write the radii of the one center projector
990!> \param atomic_kind_set ...
991!> \param qs_kind_set ...
992!> \param subsys_section ...
993!> \version 1.0
994! **************************************************************************************************
995 SUBROUTINE write_paw_radii(atomic_kind_set, qs_kind_set, subsys_section)
996 TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
997 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
998 TYPE(section_vals_type), POINTER :: subsys_section
999
1000 CHARACTER(LEN=default_string_length) :: name, unit_str
1001 INTEGER :: ikind, nkind, output_unit
1002 LOGICAL :: paw_atom
1003 REAL(kind=dp) :: conv, rcprj
1004 TYPE(cp_logger_type), POINTER :: logger
1005 TYPE(paw_proj_set_type), POINTER :: paw_proj_set
1006
1007 NULLIFY (logger)
1008 logger => cp_get_default_logger()
1009 output_unit = cp_print_key_unit_nr(logger, subsys_section, &
1010 "PRINT%RADII/GAPW_PRJ_RADII", extension=".Log")
1011 CALL section_vals_val_get(subsys_section, "PRINT%RADII%UNIT", c_val=unit_str)
1012 conv = cp_unit_from_cp2k(1.0_dp, trim(unit_str))
1013 IF (output_unit > 0) THEN
1014 nkind = SIZE(qs_kind_set)
1015 WRITE (unit=output_unit, fmt="(/,T2,A,T56,A,T63,A,T75,A)") &
1016 "RADII: ONE CENTER PROJECTORS in "// &
1017 trim(unit_str), "Kind", "Label", "Radius"
1018 DO ikind = 1, nkind
1019 CALL get_atomic_kind(atomic_kind_set(ikind), name=name)
1020 CALL get_qs_kind(qs_kind_set(ikind), &
1021 paw_atom=paw_atom, paw_proj_set=paw_proj_set)
1022 IF (paw_atom .AND. ASSOCIATED(paw_proj_set)) THEN
1023 CALL get_paw_proj_set(paw_proj_set=paw_proj_set, &
1024 rcprj=rcprj)
1025 WRITE (unit=output_unit, fmt="(T55,I5,T63,A5,T69,F12.6)") &
1026 ikind, name, rcprj*conv
1027 END IF
1028 END DO
1029 END IF
1030 CALL cp_print_key_finished_output(output_unit, logger, subsys_section, &
1031 "PRINT%RADII/GAPW_PRJ_RADII")
1032 END SUBROUTINE write_paw_radii
1033
1034END MODULE qs_interactions
All kind of helpful little routines.
Definition ao_util.F:14
real(kind=dp) function, public exp_radius(l, alpha, threshold, prefactor, epsabs, epsrel, rlow)
The radius of a primitive Gaussian function for a given threshold is calculated. g(r) = prefactor*r**...
Definition ao_util.F:96
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)
...
subroutine, public set_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, short_kind_radius, ccon)
Set the components of Gaussian-type orbital (GTO) basis set data set.
Defines control structures, which contain the parameters and the settings for the DFT-based calculati...
various routines to log and control the output. The idea is that decisions about where to log should ...
type(cp_logger_type) function, pointer, public cp_get_default_logger()
returns the default logger
routines to handle the output, The idea is to remove the decision of wheter to output and what to out...
subroutine, public cp_print_key_finished_output(unit_nr, logger, basis_section, print_key_path, local, ignore_should_output, on_file, mpi_io)
should be called after you finish working with a unit obtained with cp_print_key_unit_nr,...
integer function, public cp_print_key_unit_nr(logger, basis_section, print_key_path, extension, middle_name, local, log_filename, ignore_should_output, file_form, file_position, file_action, file_status, do_backup, on_file, is_new_file, mpi_io, fout, file_access)
...
unit conversion facility
Definition cp_units.F:30
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:1311
Definition of the atomic potential types.
objects that represent the structure of input sections and the data contained in an input section
subroutine, public section_vals_val_get(section_vals, keyword_name, i_rep_section, i_rep_val, n_rep_val, val, l_val, i_val, r_val, c_val, l_vals, i_vals, r_vals, c_vals, explicit)
returns the requested value
Defines the basic variable types.
Definition kinds.F:23
integer, parameter, public dp
Definition kinds.F:34
integer, parameter, public default_string_length
Definition kinds.F:57
Calculate MAO's and analyze wavefunctions.
Definition mao_basis.F:15
generate or use from input minimal basis set
subroutine, public get_paw_proj_set(paw_proj_set, csprj, chprj, first_prj, first_prjs, last_prj, local_oce_sphi_h, local_oce_sphi_s, maxl, ncgauprj, nsgauprj, nsatbas, nsotot, nprj, o2nindex, n2oindex, rcprj, rzetprj, zisomin, zetprj)
Get informations about a paw projectors set.
subroutine, public set_paw_proj_set(paw_proj_set, rzetprj, rcprj)
Set informations about a paw projectors set.
Types used by CNEO-DFT (see J. Chem. Theory Comput. 2025, 21, 16, 7865–7877).
Calculate the interaction radii for the operator matrix calculation.
subroutine, public write_pgf_orb_radii(basis, atomic_kind_set, qs_kind_set, subsys_section)
Write the orbital basis function radii to the output unit.
subroutine, public init_interaction_radii_orb_basis(orb_basis_set, eps_pgf_orb, eps_pgf_short)
...
subroutine, public write_paw_radii(atomic_kind_set, qs_kind_set, subsys_section)
Write the radii of the one center projector.
subroutine, public write_ppnl_radii(atomic_kind_set, qs_kind_set, subsys_section)
Write the radii of the projector functions of the Goedecker pseudopotential (GTH, non-local part) to ...
subroutine, public init_se_nlradius(se_control, atomic_kind_set, qs_kind_set, subsys_section)
...
subroutine, public write_ppl_radii(atomic_kind_set, qs_kind_set, subsys_section)
Write the radii of the exponential functions of the Goedecker pseudopotential (GTH,...
subroutine, public init_interaction_radii(qs_control, qs_kind_set)
Initialize all the atomic kind radii for a given threshold value.
subroutine, public write_core_charge_radii(atomic_kind_set, qs_kind_set, subsys_section)
Write the radii of the core charge distributions to the output unit.
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.
Utilities for string manipulations.
elemental subroutine, public uppercase(string)
Convert all lower case characters in a string to upper case.
Local UPF matrix elements from the original radial potential. The Gaussian core electrostatics suppli...
subroutine, public upf_local_init(pot, alpha, degree, local, eps)
Prepare the short-range local potential outside parallel regions.
Provides all information about an atomic kind.
type of a logger, at the moment it contains just a print level starting at which level it should be l...
Provides all information about a quickstep kind.