(git:53a46e8)
Loading...
Searching...
No Matches
qs_kind_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
8! **************************************************************************************************
9!> \brief Define the quickstep kind type and their sub types
10!> \author Ole Schuett
11!>
12!> <b>Modification history:</b>
13!> - 01.2002 creation [MK]
14!> - 04.2002 added pao [fawzi]
15!> - 09.2002 adapted for POL/KG use [GT]
16!> - 02.2004 flexible normalization of basis sets [jgh]
17!> - 03.2004 attach/detach routines [jgh]
18!> - 10.2004 removed pao [fawzi]
19!> - 08.2014 separated qs-related stuff from atomic_kind_types.F [Ole Schuett]
20!> - 07.2015 new container for basis sets [jgh]
21!> - 04.2021 init dft_plus_u_type [MK]
22! **************************************************************************************************
27 USE atom_types, ONLY: atom_ecppot_type,&
28 lmat,&
30 USE atom_upf, ONLY: atom_read_upf,&
40 USE basis_set_types, ONLY: &
51 USE cp_output_handling, ONLY: cp_p_file,&
55 USE external_potential_types, ONLY: &
60 USE input_constants, ONLY: &
69 USE kinds, ONLY: default_path_length,&
71 dp
72 USE mathconstants, ONLY: pi
75 nco,&
76 ncoset
83 ptable
84 USE physcon, ONLY: angstrom,&
85 bohr,&
86 evolt
112 USE string_utilities, ONLY: uppercase
113 USE xtb_parameters, ONLY: xtb_set_kab
118#include "./base/base_uses.f90"
119
120 IMPLICIT NONE
121
122 PRIVATE
123
124 ! Global parameters (only in this module)
125
126 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'qs_kind_types'
127
128! **************************************************************************************************
129!> \brief Input parameters for the DFT+U method
130! **************************************************************************************************
131 TYPE dft_plus_u_type
132 INTEGER :: l = -1
133 INTEGER :: n = -1
134 INTEGER :: max_scf = -1
135 REAL(KIND=dp) :: eps_u_ramping = 0.0_dp
136 REAL(KIND=dp) :: eps_scf = huge(0.0_dp)
137 REAL(KIND=dp) :: u_minus_j_target = 0.0_dp
138 REAL(KIND=dp) :: u_minus_j = 0.0_dp
139 REAL(KIND=dp) :: u_ramping = 0.0_dp
140 REAL(KIND=dp) :: u = 0.0_dp
141 REAL(KIND=dp) :: j = 0.0_dp
142 REAL(KIND=dp) :: alpha = 0.0_dp
143 REAL(KIND=dp) :: beta = 0.0_dp
144 REAL(KIND=dp) :: j0 = 0.0_dp
145 REAL(KIND=dp) :: occupation = -1.0_dp
146 INTEGER, DIMENSION(:), POINTER :: orbitals => null()
147 LOGICAL :: init_u_ramping_each_scf = .false.
148 LOGICAL :: smear = .false.
149 REAL(KIND=dp), DIMENSION(:), POINTER :: nelec => null()
150 END TYPE dft_plus_u_type
151
152! **************************************************************************************************
153!> \brief Holds information about a PAO potential
154! **************************************************************************************************
156 INTEGER :: maxl = -1
157 REAL(kind=dp) :: beta = 0.0_dp
158 REAL(kind=dp) :: weight = 0.0_dp
159 INTEGER :: max_projector = -1
160 REAL(kind=dp) :: beta_radius = huge(dp)
161 END TYPE pao_potential_type
162
163! **************************************************************************************************
164!> \brief Holds information about a PAO descriptor
165! **************************************************************************************************
167 REAL(kind=dp) :: beta = 0.0_dp
168 REAL(kind=dp) :: beta_radius = huge(dp)
169 REAL(kind=dp) :: weight = 0.0_dp
170 REAL(kind=dp) :: screening = 0.0_dp
171 REAL(kind=dp) :: screening_radius = huge(dp)
172 END TYPE pao_descriptor_type
173
174! **************************************************************************************************
175!> \brief Provides all information about a quickstep kind
176! **************************************************************************************************
178 CHARACTER(LEN=default_string_length) :: name = ""
179 CHARACTER(LEN=2) :: element_symbol = ""
180 INTEGER :: natom = -1
181 TYPE(all_potential_type), POINTER :: all_potential => null()
182 TYPE(local_potential_type), POINTER :: tnadd_potential => null()
183 TYPE(gth_potential_type), POINTER :: gth_potential => null()
184 TYPE(sgp_potential_type), POINTER :: sgp_potential => null()
185 TYPE(semi_empirical_type), POINTER :: se_parameter => null()
186 TYPE(qs_dftb_atom_type), POINTER :: dftb_parameter => null()
187 TYPE(xtb_atom_type), POINTER :: xtb_parameter => null()
188 !
189 TYPE(atom_upfpot_type), POINTER :: upf_potential => null()
190 TYPE(cneo_potential_type), POINTER :: cneo_potential => null()
191 !
193 DIMENSION(20) :: basis_sets = basis_set_container_type()
194 ! Atomic radii
195 REAL(kind=dp) :: covalent_radius = 0.0_dp
196 REAL(kind=dp) :: vdw_radius = 0.0_dp
197 ! GAPW specific data
198 TYPE(paw_proj_set_type), POINTER :: paw_proj_set => null()
199 REAL(kind=dp) :: hard_radius = 0.8_dp*bohr ! for hard and soft exp
200 REAL(kind=dp) :: hard0_radius = 0.8_dp*bohr ! for hard exp of rho0
201 REAL(kind=dp) :: max_rad_local = 13.2_dp*bohr ! max GTO radius used in GAPW
202 LOGICAL :: paw_atom = .false. ! needs atomic rho1
203 LOGICAL :: gpw_type_forced = .false. ! gpw atom even if with hard exponents
204 !
205 LOGICAL :: ghost = .false.
206 LOGICAL :: monovalent = .false.
207 LOGICAL :: floating = .false.
208 INTEGER :: lmax_dftb = -1
209 REAL(kind=dp) :: dudq_dftb3 = 0.0_dp
210 REAL(kind=dp) :: magnetization = 0.0_dp
211 INTEGER, DIMENSION(:, :), POINTER :: addel => null()
212 INTEGER, DIMENSION(:, :), POINTER :: laddel => null()
213 INTEGER, DIMENSION(:, :), POINTER :: naddel => null()
214 TYPE(harmonics_atom_type), POINTER :: harmonics => null()
215 TYPE(grid_atom_type), POINTER :: grid_atom => null()
216 INTEGER :: ngrid_rad = 50
217 INTEGER :: ngrid_ang = 50
218 INTEGER :: lmax_rho0 = 0
219 INTEGER :: mao = -1
220 INTEGER, DIMENSION(:), POINTER :: elec_conf => null() ! used to set up the initial atomic guess
221 LOGICAL :: bs_occupation = .false.
222 TYPE(dft_plus_u_type), POINTER :: dft_plus_u => null()
223 LOGICAL :: no_optimize = .true.
224 !
225 REAL(kind=dp), DIMENSION(:, :), POINTER :: nlcc_pot => null()
226 !
227 TYPE(qs_atom_dispersion_type), POINTER :: dispersion => null()
228 REAL(kind=dp), DIMENSION(:, :), POINTER :: reltmat => null()
229 INTEGER :: pao_basis_size = -1
230 CHARACTER(LEN=default_path_length) :: pao_model_file = ""
231 TYPE(pao_potential_type), DIMENSION(:), POINTER :: pao_potentials => null()
232 TYPE(pao_descriptor_type), DIMENSION(:), POINTER :: pao_descriptors => null()
233 END TYPE qs_kind_type
234
235! **************************************************************************************************
236!> \brief Provides a vector of pointers of type qs_kind_type
237! **************************************************************************************************
238 TYPE qs_kind_p_type
239 TYPE(qs_kind_type), DIMENSION(:), &
240 POINTER :: qs_kind_set => null()
241 END TYPE qs_kind_p_type
242
243 ! Public subroutines
244
245 PUBLIC :: check_qs_kind_set, &
247 get_qs_kind, &
249 has_nlcc, &
254 set_qs_kind, &
259
260 ! Public data types
262
263CONTAINS
264
265! **************************************************************************************************
266!> \brief Destructor routine for a set of qs kinds
267!> \param qs_kind_set ...
268!> \date 02.01.2002
269!> \author Matthias Krack (MK)
270!> \version 2.0
271! **************************************************************************************************
272 SUBROUTINE deallocate_qs_kind_set(qs_kind_set)
273
274 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
275
276 INTEGER :: ikind, nkind
277
278 IF (ASSOCIATED(qs_kind_set)) THEN
279
280 nkind = SIZE(qs_kind_set)
281
282 DO ikind = 1, nkind
283 IF (ASSOCIATED(qs_kind_set(ikind)%all_potential)) THEN
284 CALL deallocate_potential(qs_kind_set(ikind)%all_potential)
285 END IF
286 IF (ASSOCIATED(qs_kind_set(ikind)%tnadd_potential)) THEN
287 CALL deallocate_potential(qs_kind_set(ikind)%tnadd_potential)
288 END IF
289 IF (ASSOCIATED(qs_kind_set(ikind)%gth_potential)) THEN
290 CALL deallocate_potential(qs_kind_set(ikind)%gth_potential)
291 END IF
292 IF (ASSOCIATED(qs_kind_set(ikind)%sgp_potential)) THEN
293 CALL deallocate_potential(qs_kind_set(ikind)%sgp_potential)
294 END IF
295 IF (ASSOCIATED(qs_kind_set(ikind)%upf_potential)) THEN
296 CALL atom_release_upf(qs_kind_set(ikind)%upf_potential)
297 DEALLOCATE (qs_kind_set(ikind)%upf_potential)
298 END IF
299 IF (ASSOCIATED(qs_kind_set(ikind)%cneo_potential)) THEN
300 CALL deallocate_cneo_potential(qs_kind_set(ikind)%cneo_potential)
301 END IF
302 IF (ASSOCIATED(qs_kind_set(ikind)%se_parameter)) THEN
303 CALL semi_empirical_release(qs_kind_set(ikind)%se_parameter)
304 END IF
305 IF (ASSOCIATED(qs_kind_set(ikind)%dftb_parameter)) THEN
306 CALL deallocate_dftb_atom_param(qs_kind_set(ikind)%dftb_parameter)
307 END IF
308 IF (ASSOCIATED(qs_kind_set(ikind)%xtb_parameter)) THEN
309 CALL deallocate_xtb_atom_param(qs_kind_set(ikind)%xtb_parameter)
310 END IF
311 IF (ASSOCIATED(qs_kind_set(ikind)%paw_proj_set)) THEN
312 CALL deallocate_paw_proj_set(qs_kind_set(ikind)%paw_proj_set)
313 END IF
314 IF (ASSOCIATED(qs_kind_set(ikind)%harmonics)) THEN
315 CALL deallocate_harmonics_atom(qs_kind_set(ikind)%harmonics)
316 END IF
317 IF (ASSOCIATED(qs_kind_set(ikind)%grid_atom)) THEN
318 CALL deallocate_grid_atom(qs_kind_set(ikind)%grid_atom)
319 END IF
320 IF (ASSOCIATED(qs_kind_set(ikind)%elec_conf)) THEN
321 DEALLOCATE (qs_kind_set(ikind)%elec_conf)
322 END IF
323
324 IF (ASSOCIATED(qs_kind_set(ikind)%dft_plus_u)) THEN
325 IF (ASSOCIATED(qs_kind_set(ikind)%dft_plus_u%orbitals)) THEN
326 DEALLOCATE (qs_kind_set(ikind)%dft_plus_u%orbitals)
327 END IF
328 IF (ASSOCIATED(qs_kind_set(ikind)%dft_plus_u%nelec)) THEN
329 DEALLOCATE (qs_kind_set(ikind)%dft_plus_u%nelec)
330 END IF
331 DEALLOCATE (qs_kind_set(ikind)%dft_plus_u)
332 END IF
333
334 IF (ASSOCIATED(qs_kind_set(ikind)%nlcc_pot)) THEN
335 DEALLOCATE (qs_kind_set(ikind)%nlcc_pot)
336 END IF
337
338 IF (ASSOCIATED(qs_kind_set(ikind)%dispersion)) THEN
339 DEALLOCATE (qs_kind_set(ikind)%dispersion)
340 END IF
341 IF (ASSOCIATED(qs_kind_set(ikind)%addel)) THEN
342 DEALLOCATE (qs_kind_set(ikind)%addel)
343 END IF
344 IF (ASSOCIATED(qs_kind_set(ikind)%naddel)) THEN
345 DEALLOCATE (qs_kind_set(ikind)%naddel)
346 END IF
347 IF (ASSOCIATED(qs_kind_set(ikind)%laddel)) THEN
348 DEALLOCATE (qs_kind_set(ikind)%laddel)
349 END IF
350 IF (ASSOCIATED(qs_kind_set(ikind)%reltmat)) THEN
351 DEALLOCATE (qs_kind_set(ikind)%reltmat)
352 END IF
353
354 IF (ASSOCIATED(qs_kind_set(ikind)%pao_potentials)) THEN
355 DEALLOCATE (qs_kind_set(ikind)%pao_potentials)
356 END IF
357 IF (ASSOCIATED(qs_kind_set(ikind)%pao_descriptors)) THEN
358 DEALLOCATE (qs_kind_set(ikind)%pao_descriptors)
359 END IF
360
361 CALL remove_basis_set_container(qs_kind_set(ikind)%basis_sets)
362
363 END DO
364 DEALLOCATE (qs_kind_set)
365 ELSE
366 CALL cp_abort(__location__, &
367 "The pointer qs_kind_set is not associated and "// &
368 "cannot be deallocated")
369 END IF
370
371 END SUBROUTINE deallocate_qs_kind_set
372
373! **************************************************************************************************
374!> \brief Get attributes of an atomic kind.
375!> \param qs_kind ...
376!> \param basis_set ...
377!> \param basis_type ...
378!> \param ncgf ...
379!> \param nsgf ...
380!> \param all_potential ...
381!> \param tnadd_potential ...
382!> \param gth_potential ...
383!> \param sgp_potential ...
384!> \param upf_potential ...
385!> \param cneo_potential ...
386!> \param se_parameter ...
387!> \param dftb_parameter ...
388!> \param xtb_parameter ...
389!> \param dftb3_param ...
390!> \param zatom ...
391!> \param zeff ...
392!> \param elec_conf ...
393!> \param mao ...
394!> \param lmax_dftb ...
395!> \param alpha_core_charge ...
396!> \param ccore_charge ...
397!> \param core_charge ...
398!> \param core_charge_radius ...
399!> \param paw_proj_set ...
400!> \param paw_atom ...
401!> \param hard_radius ...
402!> \param hard0_radius ...
403!> \param max_rad_local ...
404!> \param covalent_radius ...
405!> \param vdw_radius ...
406!> \param gpw_type_forced ...
407!> \param harmonics ...
408!> \param max_iso_not0 ...
409!> \param max_s_harm ...
410!> \param grid_atom ...
411!> \param ngrid_ang ...
412!> \param ngrid_rad ...
413!> \param lmax_rho0 ...
414!> \param dft_plus_u_atom ...
415!> \param l_of_dft_plus_u ...
416!> \param n_of_dft_plus_u ...
417!> \param u_minus_j ...
418!> \param U_of_dft_plus_u ...
419!> \param J_of_dft_plus_u ...
420!> \param alpha_of_dft_plus_u ...
421!> \param beta_of_dft_plus_u ...
422!> \param J0_of_dft_plus_u ...
423!> \param occupation_of_dft_plus_u ...
424!> \param dispersion ...
425!> \param bs_occupation ...
426!> \param magnetization ...
427!> \param no_optimize ...
428!> \param addel ...
429!> \param laddel ...
430!> \param naddel ...
431!> \param orbitals ...
432!> \param max_scf ...
433!> \param eps_scf ...
434!> \param smear ...
435!> \param u_ramping ...
436!> \param u_minus_j_target ...
437!> \param eps_u_ramping ...
438!> \param init_u_ramping_each_scf ...
439!> \param reltmat ...
440!> \param ghost ...
441!> \param monovalent ...
442!> \param floating ...
443!> \param name ...
444!> \param element_symbol ...
445!> \param pao_basis_size ...
446!> \param pao_model_file ...
447!> \param pao_potentials ...
448!> \param pao_descriptors ...
449!> \param nelec ...
450! **************************************************************************************************
451 SUBROUTINE get_qs_kind(qs_kind, &
452 basis_set, basis_type, ncgf, nsgf, &
453 all_potential, tnadd_potential, gth_potential, sgp_potential, upf_potential, &
454 cneo_potential, se_parameter, dftb_parameter, xtb_parameter, &
455 dftb3_param, zatom, zeff, elec_conf, mao, lmax_dftb, &
456 alpha_core_charge, ccore_charge, core_charge, core_charge_radius, &
457 paw_proj_set, paw_atom, hard_radius, hard0_radius, max_rad_local, &
458 covalent_radius, vdw_radius, &
459 gpw_type_forced, harmonics, max_iso_not0, max_s_harm, grid_atom, &
460 ngrid_ang, ngrid_rad, lmax_rho0, &
461 dft_plus_u_atom, l_of_dft_plus_u, n_of_dft_plus_u, &
462 u_minus_j, U_of_dft_plus_u, J_of_dft_plus_u, &
463 alpha_of_dft_plus_u, beta_of_dft_plus_u, J0_of_dft_plus_u, occupation_of_dft_plus_u, dispersion, &
464 bs_occupation, magnetization, no_optimize, addel, laddel, naddel, orbitals, &
465 max_scf, eps_scf, smear, u_ramping, u_minus_j_target, eps_u_ramping, &
466 init_u_ramping_each_scf, reltmat, ghost, monovalent, floating, name, element_symbol, &
467 pao_basis_size, pao_model_file, pao_potentials, pao_descriptors, nelec)
468
469 TYPE(qs_kind_type) :: qs_kind
470 TYPE(gto_basis_set_type), OPTIONAL, POINTER :: basis_set
471 CHARACTER(len=*), OPTIONAL :: basis_type
472 INTEGER, INTENT(OUT), OPTIONAL :: ncgf, nsgf
473 TYPE(all_potential_type), OPTIONAL, POINTER :: all_potential
474 TYPE(local_potential_type), OPTIONAL, POINTER :: tnadd_potential
475 TYPE(gth_potential_type), OPTIONAL, POINTER :: gth_potential
476 TYPE(sgp_potential_type), OPTIONAL, POINTER :: sgp_potential
477 TYPE(atom_upfpot_type), OPTIONAL, POINTER :: upf_potential
478 TYPE(cneo_potential_type), OPTIONAL, POINTER :: cneo_potential
479 TYPE(semi_empirical_type), OPTIONAL, POINTER :: se_parameter
480 TYPE(qs_dftb_atom_type), OPTIONAL, POINTER :: dftb_parameter
481 TYPE(xtb_atom_type), OPTIONAL, POINTER :: xtb_parameter
482 REAL(kind=dp), INTENT(OUT), OPTIONAL :: dftb3_param
483 INTEGER, INTENT(OUT), OPTIONAL :: zatom
484 REAL(kind=dp), INTENT(OUT), OPTIONAL :: zeff
485 INTEGER, DIMENSION(:), OPTIONAL, POINTER :: elec_conf
486 INTEGER, INTENT(OUT), OPTIONAL :: mao, lmax_dftb
487 REAL(kind=dp), INTENT(OUT), OPTIONAL :: alpha_core_charge, ccore_charge, &
488 core_charge, core_charge_radius
489 TYPE(paw_proj_set_type), OPTIONAL, POINTER :: paw_proj_set
490 LOGICAL, INTENT(OUT), OPTIONAL :: paw_atom
491 REAL(kind=dp), INTENT(OUT), OPTIONAL :: hard_radius, hard0_radius, &
492 max_rad_local, covalent_radius, &
493 vdw_radius
494 LOGICAL, INTENT(OUT), OPTIONAL :: gpw_type_forced
495 TYPE(harmonics_atom_type), OPTIONAL, POINTER :: harmonics
496 INTEGER, INTENT(OUT), OPTIONAL :: max_iso_not0, max_s_harm
497 TYPE(grid_atom_type), OPTIONAL, POINTER :: grid_atom
498 INTEGER, INTENT(OUT), OPTIONAL :: ngrid_ang, ngrid_rad, lmax_rho0
499 LOGICAL, INTENT(OUT), OPTIONAL :: dft_plus_u_atom
500 INTEGER, INTENT(OUT), OPTIONAL :: l_of_dft_plus_u, n_of_dft_plus_u
501 REAL(kind=dp), INTENT(OUT), OPTIONAL :: u_minus_j, u_of_dft_plus_u, j_of_dft_plus_u, &
502 alpha_of_dft_plus_u, beta_of_dft_plus_u, j0_of_dft_plus_u, occupation_of_dft_plus_u
503 TYPE(qs_atom_dispersion_type), OPTIONAL, POINTER :: dispersion
504 LOGICAL, INTENT(OUT), OPTIONAL :: bs_occupation
505 REAL(kind=dp), INTENT(OUT), OPTIONAL :: magnetization
506 LOGICAL, INTENT(OUT), OPTIONAL :: no_optimize
507 INTEGER, DIMENSION(:, :), OPTIONAL, POINTER :: addel, laddel, naddel
508 INTEGER, DIMENSION(:), OPTIONAL, POINTER :: orbitals
509 INTEGER, OPTIONAL :: max_scf
510 REAL(kind=dp), OPTIONAL :: eps_scf
511 LOGICAL, OPTIONAL :: smear
512 REAL(kind=dp), INTENT(OUT), OPTIONAL :: u_ramping, u_minus_j_target, &
513 eps_u_ramping
514 LOGICAL, OPTIONAL :: init_u_ramping_each_scf
515 REAL(kind=dp), DIMENSION(:, :), OPTIONAL, POINTER :: reltmat
516 LOGICAL, OPTIONAL :: ghost
517 LOGICAL, INTENT(OUT), OPTIONAL :: monovalent
518 LOGICAL, OPTIONAL :: floating
519 CHARACTER(LEN=default_string_length), &
520 INTENT(OUT), OPTIONAL :: name
521 CHARACTER(LEN=2), INTENT(OUT), OPTIONAL :: element_symbol
522 INTEGER, INTENT(OUT), OPTIONAL :: pao_basis_size
523 CHARACTER(LEN=default_path_length), INTENT(OUT), &
524 OPTIONAL :: pao_model_file
525 TYPE(pao_potential_type), DIMENSION(:), OPTIONAL, &
526 POINTER :: pao_potentials
527 TYPE(pao_descriptor_type), DIMENSION(:), &
528 OPTIONAL, POINTER :: pao_descriptors
529 REAL(kind=dp), DIMENSION(:), OPTIONAL, POINTER :: nelec
530
531 CHARACTER(LEN=default_string_length) :: my_basis_type
532 INTEGER :: l
533 LOGICAL :: found
534 TYPE(gto_basis_set_type), POINTER :: tmp_basis_set
535
536 ! Retrieve basis set from the kind container
537 IF (PRESENT(basis_type)) THEN
538 my_basis_type = basis_type
539 ELSE
540 my_basis_type = "ORB"
541 END IF
542
543 IF (PRESENT(basis_set)) THEN
544 CALL get_basis_from_container(qs_kind%basis_sets, basis_set=basis_set, &
545 basis_type=my_basis_type)
546 END IF
547
548 IF (PRESENT(ncgf)) THEN
549 CALL get_basis_from_container(qs_kind%basis_sets, basis_set=tmp_basis_set, &
550 basis_type=my_basis_type)
551 IF (ASSOCIATED(tmp_basis_set)) THEN
552 CALL get_gto_basis_set(gto_basis_set=tmp_basis_set, ncgf=ncgf)
553 ELSE IF (ASSOCIATED(qs_kind%dftb_parameter)) THEN
554 l = qs_kind%dftb_parameter%lmax
555 ncgf = ((l + 1)*(l + 2)*(l + 3))/6
556 ELSE
557 ncgf = 0
558 END IF
559 END IF
560
561 IF (PRESENT(nsgf)) THEN
562 CALL get_basis_from_container(qs_kind%basis_sets, basis_set=tmp_basis_set, &
563 basis_type=my_basis_type)
564 IF (ASSOCIATED(tmp_basis_set)) THEN
565 CALL get_gto_basis_set(gto_basis_set=tmp_basis_set, nsgf=nsgf)
566 ELSE IF (ASSOCIATED(qs_kind%dftb_parameter)) THEN
567 nsgf = qs_kind%dftb_parameter%natorb
568 ELSE
569 nsgf = 0
570 END IF
571 END IF
572
573 IF (PRESENT(all_potential)) all_potential => qs_kind%all_potential
574 IF (PRESENT(tnadd_potential)) tnadd_potential => qs_kind%tnadd_potential
575 IF (PRESENT(gth_potential)) gth_potential => qs_kind%gth_potential
576 IF (PRESENT(sgp_potential)) sgp_potential => qs_kind%sgp_potential
577 IF (PRESENT(upf_potential)) upf_potential => qs_kind%upf_potential
578 IF (PRESENT(cneo_potential)) cneo_potential => qs_kind%cneo_potential
579 IF (PRESENT(se_parameter)) se_parameter => qs_kind%se_parameter
580 IF (PRESENT(dftb_parameter)) dftb_parameter => qs_kind%dftb_parameter
581 IF (PRESENT(xtb_parameter)) xtb_parameter => qs_kind%xtb_parameter
582 IF (PRESENT(element_symbol)) element_symbol = qs_kind%element_symbol
583 IF (PRESENT(name)) name = qs_kind%name
584 IF (PRESENT(dftb3_param)) dftb3_param = qs_kind%dudq_dftb3
585 IF (PRESENT(elec_conf)) elec_conf => qs_kind%elec_conf
586 IF (PRESENT(alpha_core_charge)) THEN
587 IF (ASSOCIATED(qs_kind%all_potential)) THEN
588 CALL get_potential(potential=qs_kind%all_potential, &
589 alpha_core_charge=alpha_core_charge)
590 ELSE IF (ASSOCIATED(qs_kind%gth_potential)) THEN
591 CALL get_potential(potential=qs_kind%gth_potential, &
592 alpha_core_charge=alpha_core_charge)
593 ELSE IF (ASSOCIATED(qs_kind%sgp_potential)) THEN
594 CALL get_potential(potential=qs_kind%sgp_potential, &
595 alpha_core_charge=alpha_core_charge)
596 ELSE IF (ASSOCIATED(qs_kind%cneo_potential)) THEN
597 cpabort("CNEO ALPHA CORE CHARGE NOT AVAILABLE")
598 ELSE
599 alpha_core_charge = 1.0_dp
600 END IF
601 END IF
602 IF (PRESENT(ccore_charge)) THEN
603 IF (ASSOCIATED(qs_kind%all_potential)) THEN
604 CALL get_potential(potential=qs_kind%all_potential, &
605 ccore_charge=ccore_charge)
606 ELSE IF (ASSOCIATED(qs_kind%gth_potential)) THEN
607 CALL get_potential(potential=qs_kind%gth_potential, &
608 ccore_charge=ccore_charge)
609 ELSE IF (ASSOCIATED(qs_kind%sgp_potential)) THEN
610 CALL get_potential(potential=qs_kind%sgp_potential, &
611 ccore_charge=ccore_charge)
612 ELSE IF (ASSOCIATED(qs_kind%upf_potential)) THEN
613 cpabort("UPF CCORE CHARGE NOT AVAILABLE")
614 ELSE IF (ASSOCIATED(qs_kind%cneo_potential)) THEN
615 cpabort("CNEO CCORE CHARGE NOT AVAILABLE")
616 ELSE
617 ccore_charge = 0.0_dp
618 END IF
619 END IF
620 IF (PRESENT(core_charge_radius)) THEN
621 IF (ASSOCIATED(qs_kind%all_potential)) THEN
622 CALL get_potential(potential=qs_kind%all_potential, &
623 core_charge_radius=core_charge_radius)
624 ELSE IF (ASSOCIATED(qs_kind%gth_potential)) THEN
625 CALL get_potential(potential=qs_kind%gth_potential, &
626 core_charge_radius=core_charge_radius)
627 ELSE IF (ASSOCIATED(qs_kind%sgp_potential)) THEN
628 CALL get_potential(potential=qs_kind%sgp_potential, &
629 core_charge_radius=core_charge_radius)
630 ELSE IF (ASSOCIATED(qs_kind%upf_potential)) THEN
631 cpabort("UPF CORE CHARGE RADIUS NOT AVAILABLE")
632 ELSE IF (ASSOCIATED(qs_kind%cneo_potential)) THEN
633 cpabort("CNEO CORE CHARGE RADIUS NOT AVAILABLE")
634 ELSE
635 core_charge_radius = 0.0_dp
636 END IF
637 END IF
638 IF (PRESENT(core_charge)) THEN
639 IF (ASSOCIATED(qs_kind%all_potential)) THEN
640 CALL get_potential(potential=qs_kind%all_potential, &
641 zeff=core_charge)
642 ELSE IF (ASSOCIATED(qs_kind%gth_potential)) THEN
643 CALL get_potential(potential=qs_kind%gth_potential, &
644 zeff=core_charge)
645 ELSE IF (ASSOCIATED(qs_kind%sgp_potential)) THEN
646 CALL get_potential(potential=qs_kind%sgp_potential, &
647 zeff=core_charge)
648 ELSE IF (ASSOCIATED(qs_kind%upf_potential)) THEN
649 cpabort("UPF CORE CHARGE NOT AVAILABLE")
650 ELSE IF (ASSOCIATED(qs_kind%cneo_potential)) THEN
651 CALL get_cneo_potential(potential=qs_kind%cneo_potential, &
652 zeff=core_charge)
653 ELSE
654 core_charge = 0.0_dp
655 END IF
656 END IF
657
658 IF (PRESENT(zatom)) THEN
659 ! Retrieve information on element
660 CALL get_ptable_info(qs_kind%element_symbol, ielement=zatom, found=found)
661 cpassert(found)
662 END IF
663
664 IF (PRESENT(zeff)) THEN
665 IF (ASSOCIATED(qs_kind%all_potential)) THEN
666 CALL get_potential(potential=qs_kind%all_potential, zeff=zeff)
667 ELSE IF (ASSOCIATED(qs_kind%gth_potential)) THEN
668 CALL get_potential(potential=qs_kind%gth_potential, zeff=zeff)
669 ELSE IF (ASSOCIATED(qs_kind%sgp_potential)) THEN
670 CALL get_potential(potential=qs_kind%sgp_potential, zeff=zeff)
671 ELSE IF (ASSOCIATED(qs_kind%upf_potential)) THEN
672 zeff = qs_kind%upf_potential%zion
673 ELSE IF (ASSOCIATED(qs_kind%cneo_potential)) THEN
674 CALL get_cneo_potential(potential=qs_kind%cneo_potential, zeff=zeff)
675 ELSE
676 zeff = 0.0_dp
677 END IF
678 END IF
679
680 IF (PRESENT(covalent_radius)) covalent_radius = qs_kind%covalent_radius
681 IF (PRESENT(vdw_radius)) vdw_radius = qs_kind%vdw_radius
682
683 IF (PRESENT(paw_proj_set)) paw_proj_set => qs_kind%paw_proj_set
684 IF (PRESENT(paw_atom)) paw_atom = qs_kind%paw_atom
685 IF (PRESENT(gpw_type_forced)) gpw_type_forced = qs_kind%gpw_type_forced
686 IF (PRESENT(hard_radius)) hard_radius = qs_kind%hard_radius
687 IF (PRESENT(hard0_radius)) hard0_radius = qs_kind%hard0_radius
688 IF (PRESENT(max_rad_local)) max_rad_local = qs_kind%max_rad_local
689 IF (PRESENT(harmonics)) harmonics => qs_kind%harmonics
690 IF (PRESENT(max_s_harm)) THEN
691 IF (ASSOCIATED(qs_kind%harmonics)) THEN
692 max_s_harm = qs_kind%harmonics%max_s_harm
693 ELSE
694 max_s_harm = 0
695 END IF
696 END IF
697 IF (PRESENT(max_iso_not0)) THEN
698 IF (ASSOCIATED(qs_kind%harmonics)) THEN
699 max_iso_not0 = qs_kind%harmonics%max_iso_not0
700 ELSE
701 max_iso_not0 = 0
702 END IF
703 END IF
704 IF (PRESENT(grid_atom)) grid_atom => qs_kind%grid_atom
705 IF (PRESENT(ngrid_ang)) ngrid_ang = qs_kind%ngrid_ang
706 IF (PRESENT(ngrid_rad)) ngrid_rad = qs_kind%ngrid_rad
707 IF (PRESENT(lmax_rho0)) lmax_rho0 = qs_kind%lmax_rho0
708 IF (PRESENT(ghost)) ghost = qs_kind%ghost
709 IF (PRESENT(monovalent)) monovalent = qs_kind%monovalent
710 IF (PRESENT(floating)) floating = qs_kind%floating
711 IF (PRESENT(dft_plus_u_atom)) dft_plus_u_atom = ASSOCIATED(qs_kind%dft_plus_u)
712 IF (PRESENT(l_of_dft_plus_u)) THEN
713 IF (ASSOCIATED(qs_kind%dft_plus_u)) THEN
714 l_of_dft_plus_u = qs_kind%dft_plus_u%l
715 ELSE
716 l_of_dft_plus_u = -1
717 END IF
718 END IF
719 IF (PRESENT(n_of_dft_plus_u)) THEN
720 IF (ASSOCIATED(qs_kind%dft_plus_u)) THEN
721 n_of_dft_plus_u = qs_kind%dft_plus_u%n
722 ELSE
723 n_of_dft_plus_u = -1
724 END IF
725 END IF
726 IF (PRESENT(u_minus_j)) THEN
727 IF (ASSOCIATED(qs_kind%dft_plus_u)) THEN
728 u_minus_j = qs_kind%dft_plus_u%u_minus_j
729 ELSE
730 u_minus_j = 0.0_dp
731 END IF
732 END IF
733 IF (PRESENT(u_minus_j_target)) THEN
734 IF (ASSOCIATED(qs_kind%dft_plus_u)) THEN
735 u_minus_j_target = qs_kind%dft_plus_u%u_minus_j_target
736 ELSE
737 u_minus_j_target = 0.0_dp
738 END IF
739 END IF
740 IF (PRESENT(u_of_dft_plus_u)) THEN
741 IF (ASSOCIATED(qs_kind%dft_plus_u)) THEN
742 u_of_dft_plus_u = qs_kind%dft_plus_u%U
743 ELSE
744 u_of_dft_plus_u = 0.0_dp
745 END IF
746 END IF
747 IF (PRESENT(j_of_dft_plus_u)) THEN
748 IF (ASSOCIATED(qs_kind%dft_plus_u)) THEN
749 j_of_dft_plus_u = qs_kind%dft_plus_u%J
750 ELSE
751 j_of_dft_plus_u = 0.0_dp
752 END IF
753 END IF
754 IF (PRESENT(alpha_of_dft_plus_u)) THEN
755 IF (ASSOCIATED(qs_kind%dft_plus_u)) THEN
756 alpha_of_dft_plus_u = qs_kind%dft_plus_u%alpha
757 ELSE
758 alpha_of_dft_plus_u = 0.0_dp
759 END IF
760 END IF
761 IF (PRESENT(beta_of_dft_plus_u)) THEN
762 IF (ASSOCIATED(qs_kind%dft_plus_u)) THEN
763 beta_of_dft_plus_u = qs_kind%dft_plus_u%beta
764 ELSE
765 beta_of_dft_plus_u = 0.0_dp
766 END IF
767 END IF
768 IF (PRESENT(j0_of_dft_plus_u)) THEN
769 IF (ASSOCIATED(qs_kind%dft_plus_u)) THEN
770 j0_of_dft_plus_u = qs_kind%dft_plus_u%J0
771 ELSE
772 j0_of_dft_plus_u = 0.0_dp
773 END IF
774 END IF
775 IF (PRESENT(occupation_of_dft_plus_u)) THEN
776 IF (ASSOCIATED(qs_kind%dft_plus_u)) THEN
777 occupation_of_dft_plus_u = qs_kind%dft_plus_u%occupation
778 ELSE
779 occupation_of_dft_plus_u = -1.0_dp
780 END IF
781 END IF
782
783 IF (PRESENT(init_u_ramping_each_scf)) THEN
784 IF (ASSOCIATED(qs_kind%dft_plus_u)) THEN
785 init_u_ramping_each_scf = qs_kind%dft_plus_u%init_u_ramping_each_scf
786 ELSE
787 init_u_ramping_each_scf = .false.
788 END IF
789 END IF
790 IF (PRESENT(u_ramping)) THEN
791 IF (ASSOCIATED(qs_kind%dft_plus_u)) THEN
792 u_ramping = qs_kind%dft_plus_u%u_ramping
793 ELSE
794 u_ramping = 0.0_dp
795 END IF
796 END IF
797 IF (PRESENT(eps_u_ramping)) THEN
798 IF (ASSOCIATED(qs_kind%dft_plus_u)) THEN
799 eps_u_ramping = qs_kind%dft_plus_u%eps_u_ramping
800 ELSE
801 eps_u_ramping = 1.0e-5_dp
802 END IF
803 END IF
804 IF (PRESENT(nelec)) THEN
805 NULLIFY (nelec)
806 IF (ASSOCIATED(qs_kind%dft_plus_u)) THEN
807 IF (ASSOCIATED(qs_kind%dft_plus_u%nelec)) THEN
808 nelec => qs_kind%dft_plus_u%nelec
809 END IF
810 END IF
811 END IF
812 IF (PRESENT(orbitals)) THEN
813 NULLIFY (orbitals)
814 IF (ASSOCIATED(qs_kind%dft_plus_u)) THEN
815 IF (ASSOCIATED(qs_kind%dft_plus_u%orbitals)) THEN
816 orbitals => qs_kind%dft_plus_u%orbitals
817 END IF
818 END IF
819 END IF
820 IF (PRESENT(eps_scf)) THEN
821 IF (ASSOCIATED(qs_kind%dft_plus_u)) THEN
822 eps_scf = qs_kind%dft_plus_u%eps_scf
823 ELSE
824 eps_scf = 1.0e30_dp
825 END IF
826 END IF
827 IF (PRESENT(max_scf)) THEN
828 IF (ASSOCIATED(qs_kind%dft_plus_u)) THEN
829 max_scf = qs_kind%dft_plus_u%max_scf
830 ELSE
831 max_scf = -1
832 END IF
833 END IF
834 IF (PRESENT(smear)) THEN
835 IF (ASSOCIATED(qs_kind%dft_plus_u)) THEN
836 smear = qs_kind%dft_plus_u%smear
837 ELSE
838 smear = .false.
839 END IF
840 END IF
841 IF (PRESENT(dispersion)) dispersion => qs_kind%dispersion
842 IF (PRESENT(bs_occupation)) bs_occupation = qs_kind%bs_occupation
843 IF (PRESENT(addel)) addel => qs_kind%addel
844 IF (PRESENT(laddel)) laddel => qs_kind%laddel
845 IF (PRESENT(naddel)) naddel => qs_kind%naddel
846
847 IF (PRESENT(magnetization)) magnetization = qs_kind%magnetization
848
849 IF (PRESENT(no_optimize)) no_optimize = qs_kind%no_optimize
850
851 IF (PRESENT(reltmat)) reltmat => qs_kind%reltmat
852
853 IF (PRESENT(mao)) mao = qs_kind%mao
854
855 IF (PRESENT(lmax_dftb)) lmax_dftb = qs_kind%lmax_dftb
856
857 IF (PRESENT(pao_basis_size)) pao_basis_size = qs_kind%pao_basis_size
858 IF (PRESENT(pao_model_file)) pao_model_file = qs_kind%pao_model_file
859 IF (PRESENT(pao_potentials)) pao_potentials => qs_kind%pao_potentials
860 IF (PRESENT(pao_descriptors)) pao_descriptors => qs_kind%pao_descriptors
861 END SUBROUTINE get_qs_kind
862
863! **************************************************************************************************
864!> \brief Get attributes of an atomic kind set.
865!> \param qs_kind_set ...
866!> \param all_potential_present ...
867!> \param tnadd_potential_present ...
868!> \param gth_potential_present ...
869!> \param sgp_potential_present ...
870!> \param paw_atom_present ...
871!> \param dft_plus_u_atom_present ...
872!> \param maxcgf ...
873!> \param maxsgf ...
874!> \param maxco ...
875!> \param maxco_proj ...
876!> \param maxgtops ...
877!> \param maxlgto ...
878!> \param maxlprj ...
879!> \param maxnset ...
880!> \param maxsgf_set ...
881!> \param ncgf ...
882!> \param npgf ...
883!> \param nset ...
884!> \param nsgf ...
885!> \param nshell ...
886!> \param maxpol ...
887!> \param maxlppl ...
888!> \param maxlppnl ...
889!> \param maxppnl ...
890!> \param nelectron ...
891!> \param maxder ...
892!> \param max_ngrid_rad ...
893!> \param max_sph_harm ...
894!> \param maxg_iso_not0 ...
895!> \param lmax_rho0 ...
896!> \param basis_rcut ...
897!> \param basis_type ...
898!> \param total_zeff_corr ... [SGh]
899!> \param npgf_seg total number of primitive GTOs in "segmented contraction format"
900!> \param cneo_potential_present ...
901!> \param nkind_q ...
902!> \param natom_q ...
903! **************************************************************************************************
904 SUBROUTINE get_qs_kind_set(qs_kind_set, &
905 all_potential_present, tnadd_potential_present, gth_potential_present, &
906 sgp_potential_present, paw_atom_present, dft_plus_u_atom_present, &
907 maxcgf, maxsgf, maxco, maxco_proj, maxgtops, maxlgto, maxlprj, maxnset, maxsgf_set, &
908 ncgf, npgf, nset, nsgf, nshell, maxpol, maxlppl, maxlppnl, maxppnl, &
909 nelectron, maxder, max_ngrid_rad, max_sph_harm, maxg_iso_not0, lmax_rho0, &
910 basis_rcut, &
911 basis_type, total_zeff_corr, npgf_seg, &
912 cneo_potential_present, nkind_q, natom_q)
913
914 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
915 LOGICAL, INTENT(OUT), OPTIONAL :: all_potential_present, tnadd_potential_present, &
916 gth_potential_present, sgp_potential_present, paw_atom_present, dft_plus_u_atom_present
917 INTEGER, INTENT(OUT), OPTIONAL :: maxcgf, maxsgf, maxco, maxco_proj, maxgtops, maxlgto, &
918 maxlprj, maxnset, maxsgf_set, ncgf, npgf, nset, nsgf, nshell, maxpol, maxlppl, maxlppnl, &
919 maxppnl, nelectron
920 INTEGER, INTENT(IN), OPTIONAL :: maxder
921 INTEGER, INTENT(OUT), OPTIONAL :: max_ngrid_rad, max_sph_harm, &
922 maxg_iso_not0, lmax_rho0
923 REAL(kind=dp), INTENT(OUT), OPTIONAL :: basis_rcut
924 CHARACTER(len=*), OPTIONAL :: basis_type
925 REAL(kind=dp), INTENT(OUT), OPTIONAL :: total_zeff_corr
926 INTEGER, INTENT(OUT), OPTIONAL :: npgf_seg
927 LOGICAL, INTENT(OUT), OPTIONAL :: cneo_potential_present
928 INTEGER, INTENT(OUT), OPTIONAL :: nkind_q, natom_q
929
930 CHARACTER(len=default_string_length) :: my_basis_type
931 INTEGER :: ikind, imax, lmax_rho0_kind, &
932 max_iso_not0, max_s_harm, n, &
933 ngrid_rad, nkind, nrloc(10), &
934 nrpot(1:15, 0:10)
935 LOGICAL :: dft_plus_u_atom, ecp_semi_local, paw_atom
936 REAL(kind=dp) :: brcut, zeff, zeff_correction
937 TYPE(all_potential_type), POINTER :: all_potential
938 TYPE(cneo_potential_type), POINTER :: cneo_potential
939 TYPE(gth_potential_type), POINTER :: gth_potential
940 TYPE(gto_basis_set_type), POINTER :: tmp_basis_set
941 TYPE(local_potential_type), POINTER :: tnadd_potential
942 TYPE(paw_proj_set_type), POINTER :: paw_proj_set
943 TYPE(qs_dftb_atom_type), POINTER :: dftb_parameter
944 TYPE(qs_kind_type), POINTER :: qs_kind
945 TYPE(sgp_potential_type), POINTER :: sgp_potential
946
947 IF (PRESENT(basis_type)) THEN
948 my_basis_type = basis_type
949 ELSE
950 my_basis_type = "ORB"
951 END IF
952
953 IF (ASSOCIATED(qs_kind_set)) THEN
954
955 IF (PRESENT(maxcgf)) maxcgf = 0
956 IF (PRESENT(maxco)) maxco = 0
957 IF (PRESENT(maxco_proj)) maxco_proj = 0
958 IF (PRESENT(maxg_iso_not0)) maxg_iso_not0 = 0
959 IF (PRESENT(maxgtops)) maxgtops = 0
960 IF (PRESENT(maxlgto)) maxlgto = -1
961 IF (PRESENT(maxlppl)) maxlppl = -1
962 IF (PRESENT(maxlppnl)) maxlppnl = -1
963 IF (PRESENT(maxpol)) maxpol = -1
964 IF (PRESENT(maxlprj)) maxlprj = -1
965 IF (PRESENT(maxnset)) maxnset = 0
966 IF (PRESENT(maxppnl)) maxppnl = 0
967 IF (PRESENT(maxsgf)) maxsgf = 0
968 IF (PRESENT(maxsgf_set)) maxsgf_set = 0
969 IF (PRESENT(ncgf)) ncgf = 0
970 IF (PRESENT(nelectron)) nelectron = 0
971 IF (PRESENT(npgf)) npgf = 0
972 IF (PRESENT(nset)) nset = 0
973 IF (PRESENT(nsgf)) nsgf = 0
974 IF (PRESENT(nshell)) nshell = 0
975 IF (PRESENT(all_potential_present)) all_potential_present = .false.
976 IF (PRESENT(tnadd_potential_present)) tnadd_potential_present = .false.
977 IF (PRESENT(gth_potential_present)) gth_potential_present = .false.
978 IF (PRESENT(sgp_potential_present)) sgp_potential_present = .false.
979 IF (PRESENT(cneo_potential_present)) cneo_potential_present = .false.
980 IF (PRESENT(nkind_q)) nkind_q = 0
981 IF (PRESENT(natom_q)) natom_q = 0
982 IF (PRESENT(paw_atom_present)) paw_atom_present = .false.
983 IF (PRESENT(max_ngrid_rad)) max_ngrid_rad = 0
984 IF (PRESENT(max_sph_harm)) max_sph_harm = 0
985 IF (PRESENT(lmax_rho0)) lmax_rho0 = 0
986 IF (PRESENT(basis_rcut)) basis_rcut = 0.0_dp
987 IF (PRESENT(total_zeff_corr)) total_zeff_corr = 0.0_dp
988 IF (PRESENT(npgf_seg)) npgf_seg = 0
989
990 nkind = SIZE(qs_kind_set)
991 DO ikind = 1, nkind
992 qs_kind => qs_kind_set(ikind)
993 CALL get_qs_kind(qs_kind=qs_kind, &
994 all_potential=all_potential, &
995 tnadd_potential=tnadd_potential, &
996 gth_potential=gth_potential, &
997 sgp_potential=sgp_potential, &
998 cneo_potential=cneo_potential, &
999 paw_proj_set=paw_proj_set, &
1000 dftb_parameter=dftb_parameter, &
1001 ngrid_rad=ngrid_rad, &
1002 max_s_harm=max_s_harm, &
1003 max_iso_not0=max_iso_not0, &
1004 paw_atom=paw_atom, &
1005 dft_plus_u_atom=dft_plus_u_atom, &
1006 lmax_rho0=lmax_rho0_kind)
1007
1008 IF (PRESENT(maxlppl) .AND. ASSOCIATED(gth_potential)) THEN
1009 CALL get_potential(potential=gth_potential, nexp_ppl=n)
1010 maxlppl = max(maxlppl, 2*(n - 1))
1011 ELSEIF (PRESENT(maxlppl) .AND. ASSOCIATED(sgp_potential)) THEN
1012 CALL get_potential(potential=sgp_potential, nrloc=nrloc, ecp_semi_local=ecp_semi_local)
1013 n = maxval(nrloc) - 2
1014 maxlppl = max(maxlppl, 2*(n - 1))
1015 IF (ecp_semi_local) THEN
1016 CALL get_potential(potential=sgp_potential, sl_lmax=imax, nrpot=nrpot)
1017 n = maxval(nrpot) - 2
1018 n = 2*(n - 1) + imax
1019 maxlppl = max(maxlppl, n)
1020 END IF
1021 END IF
1022
1023 IF (PRESENT(maxlppnl) .AND. ASSOCIATED(gth_potential)) THEN
1024 CALL get_potential(potential=gth_potential, lprj_ppnl_max=imax)
1025 maxlppnl = max(maxlppnl, imax)
1026 ELSEIF (PRESENT(maxlppnl) .AND. ASSOCIATED(sgp_potential)) THEN
1027 CALL get_potential(potential=sgp_potential, lmax=imax)
1028 maxlppnl = max(maxlppnl, imax)
1029 END IF
1030
1031 IF (PRESENT(maxpol) .AND. ASSOCIATED(tnadd_potential)) THEN
1032 CALL get_potential(potential=tnadd_potential, npol=n)
1033 maxpol = max(maxpol, 2*(n - 1))
1034 END IF
1035
1036 IF (PRESENT(maxco_proj) .AND. ASSOCIATED(paw_proj_set)) THEN
1037 CALL get_paw_proj_set(paw_proj_set=paw_proj_set, ncgauprj=imax)
1038 maxco_proj = max(maxco_proj, imax)
1039 END IF
1040
1041 IF (PRESENT(maxlprj) .AND. ASSOCIATED(paw_proj_set)) THEN
1042 CALL get_paw_proj_set(paw_proj_set=paw_proj_set, maxl=imax)
1043 maxlprj = max(maxlprj, imax)
1044 END IF
1045
1046 IF (PRESENT(maxppnl) .AND. ASSOCIATED(gth_potential)) THEN
1047 CALL get_potential(potential=gth_potential, nppnl=imax)
1048 maxppnl = max(maxppnl, imax)
1049 ELSEIF (PRESENT(maxppnl) .AND. ASSOCIATED(sgp_potential)) THEN
1050 CALL get_potential(potential=sgp_potential, nppnl=imax)
1051 maxppnl = max(maxppnl, imax)
1052 END IF
1053
1054 IF (my_basis_type(1:3) == "NUC" .AND. .NOT. ASSOCIATED(cneo_potential)) THEN
1055 NULLIFY (tmp_basis_set)
1056 ELSE
1057 CALL get_basis_from_container(qs_kind%basis_sets, basis_set=tmp_basis_set, &
1058 basis_type=my_basis_type)
1059 END IF
1060
1061 IF (PRESENT(maxcgf)) THEN
1062 IF (ASSOCIATED(tmp_basis_set)) THEN
1063 CALL get_gto_basis_set(gto_basis_set=tmp_basis_set, ncgf=imax)
1064 maxcgf = max(maxcgf, imax)
1065 ELSE IF (ASSOCIATED(qs_kind%dftb_parameter)) THEN
1066 CALL get_dftb_atom_param(dftb_parameter=dftb_parameter, lmax=imax)
1067 imax = ((imax + 1)*(imax + 2)*(imax + 3))/6
1068 maxcgf = max(maxcgf, imax)
1069 END IF
1070 END IF
1071
1072 IF (PRESENT(maxco)) THEN
1073 IF (ASSOCIATED(tmp_basis_set)) THEN
1074 IF (PRESENT(maxder)) THEN
1075 CALL get_gto_basis_set(gto_basis_set=tmp_basis_set, &
1076 maxco=imax, maxder=maxder)
1077 ELSE
1078 CALL get_gto_basis_set(gto_basis_set=tmp_basis_set, maxco=imax)
1079 END IF
1080 maxco = max(maxco, imax)
1081 END IF
1082 IF (ASSOCIATED(gth_potential)) THEN
1083 CALL get_potential(potential=gth_potential, lprj_ppnl_max=imax)
1084 maxco = max(maxco, ncoset(imax))
1085 END IF
1086 IF (ASSOCIATED(sgp_potential)) THEN
1087 CALL get_potential(potential=sgp_potential, lmax=imax)
1088 maxco = max(maxco, ncoset(imax))
1089 CALL get_potential(potential=sgp_potential, sl_lmax=imax)
1090 maxco = max(maxco, ncoset(imax))
1091 END IF
1092 END IF
1093
1094 IF (PRESENT(maxgtops)) THEN
1095 IF (ASSOCIATED(tmp_basis_set)) THEN
1096 CALL get_gto_basis_set(gto_basis_set=tmp_basis_set, maxso=imax, nset=n)
1097 maxgtops = max(maxgtops, n*imax)
1098 END IF
1099 END IF
1100
1101 IF (PRESENT(maxlgto)) THEN
1102 IF (ASSOCIATED(tmp_basis_set)) THEN
1103 CALL get_gto_basis_set(gto_basis_set=tmp_basis_set, maxl=imax)
1104 maxlgto = max(maxlgto, imax)
1105 ELSE IF (ASSOCIATED(qs_kind%dftb_parameter)) THEN
1106 CALL get_dftb_atom_param(dftb_parameter=dftb_parameter, lmax=imax)
1107 maxlgto = max(maxlgto, imax)
1108 END IF
1109 END IF
1110
1111 IF (PRESENT(maxnset)) THEN
1112 IF (ASSOCIATED(tmp_basis_set)) THEN
1113 CALL get_gto_basis_set(gto_basis_set=tmp_basis_set, nset=n)
1114 maxnset = max(maxnset, n)
1115 END IF
1116 END IF
1117
1118 IF (PRESENT(maxsgf)) THEN
1119 IF (ASSOCIATED(tmp_basis_set)) THEN
1120 CALL get_gto_basis_set(gto_basis_set=tmp_basis_set, nsgf=imax)
1121 maxsgf = max(maxsgf, imax)
1122 END IF
1123 END IF
1124
1125 IF (PRESENT(maxsgf_set)) THEN
1126 IF (ASSOCIATED(tmp_basis_set)) THEN
1127 CALL get_gto_basis_set(gto_basis_set=tmp_basis_set, maxsgf_set=imax)
1128 maxsgf_set = max(maxsgf_set, imax)
1129 END IF
1130 END IF
1131
1132 IF (PRESENT(ncgf)) THEN
1133 IF (ASSOCIATED(tmp_basis_set)) THEN
1134 CALL get_gto_basis_set(gto_basis_set=tmp_basis_set, ncgf=n)
1135 ncgf = ncgf + n*qs_kind_set(ikind)%natom
1136 ELSE IF (ASSOCIATED(qs_kind%dftb_parameter)) THEN
1137 CALL get_dftb_atom_param(dftb_parameter=dftb_parameter, lmax=imax)
1138 n = ((imax + 1)*(imax + 2)*(imax + 3))/6
1139 ncgf = ncgf + n*qs_kind_set(ikind)%natom
1140 END IF
1141 END IF
1142
1143 IF (PRESENT(npgf)) THEN
1144 IF (ASSOCIATED(tmp_basis_set)) THEN
1145 CALL get_gto_basis_set(gto_basis_set=tmp_basis_set, npgf_sum=n)
1146 npgf = npgf + n*qs_kind_set(ikind)%natom
1147 END IF
1148 END IF
1149
1150 IF (PRESENT(nset)) THEN
1151 IF (ASSOCIATED(tmp_basis_set)) THEN
1152 CALL get_gto_basis_set(gto_basis_set=tmp_basis_set, nset=n)
1153 nset = nset + n*qs_kind_set(ikind)%natom
1154 END IF
1155 END IF
1156
1157 IF (PRESENT(nsgf)) THEN
1158 IF (ASSOCIATED(tmp_basis_set)) THEN
1159 CALL get_gto_basis_set(gto_basis_set=tmp_basis_set, nsgf=n)
1160 nsgf = nsgf + n*qs_kind_set(ikind)%natom
1161 ELSE IF (ASSOCIATED(qs_kind%dftb_parameter)) THEN
1162 CALL get_dftb_atom_param(dftb_parameter=dftb_parameter, natorb=n)
1163 nsgf = nsgf + n*qs_kind_set(ikind)%natom
1164 END IF
1165 END IF
1166
1167 IF (PRESENT(nshell)) THEN
1168 IF (ASSOCIATED(tmp_basis_set)) THEN
1169 CALL get_gto_basis_set(gto_basis_set=tmp_basis_set, nshell_sum=n)
1170 nshell = nshell + n*qs_kind_set(ikind)%natom
1171 ELSE IF (ASSOCIATED(qs_kind%dftb_parameter)) THEN
1172 CALL get_dftb_atom_param(dftb_parameter=dftb_parameter, lmax=n)
1173 nshell = nshell + (n + 1)*qs_kind_set(ikind)%natom
1174 END IF
1175 END IF
1176
1177 IF (PRESENT(nelectron)) THEN
1178 IF (ASSOCIATED(qs_kind%all_potential)) THEN
1179 CALL get_potential(potential=qs_kind%all_potential, &
1180 zeff=zeff, zeff_correction=zeff_correction)
1181 ELSE IF (ASSOCIATED(qs_kind%gth_potential)) THEN
1182 CALL get_potential(potential=qs_kind%gth_potential, &
1183 zeff=zeff, zeff_correction=zeff_correction)
1184 ELSE IF (ASSOCIATED(qs_kind%sgp_potential)) THEN
1185 CALL get_potential(potential=qs_kind%sgp_potential, &
1186 zeff=zeff, zeff_correction=zeff_correction)
1187 ELSE IF (ASSOCIATED(qs_kind%cneo_potential)) THEN
1188 CALL get_cneo_potential(potential=qs_kind%cneo_potential, &
1189 zeff=zeff)
1190 zeff_correction = 0.0_dp
1191 ELSE
1192 zeff = 0.0_dp
1193 zeff_correction = 0.0_dp
1194 END IF
1195 nelectron = nelectron + qs_kind_set(ikind)%natom*nint(zeff - zeff_correction)
1196 END IF
1197
1198 IF (PRESENT(basis_rcut)) THEN
1199 IF (ASSOCIATED(tmp_basis_set)) THEN
1200 CALL get_gto_basis_set(gto_basis_set=tmp_basis_set, kind_radius=brcut)
1201 basis_rcut = max(basis_rcut, brcut)
1202 ELSE IF (ASSOCIATED(qs_kind%dftb_parameter)) THEN
1203 CALL get_dftb_atom_param(dftb_parameter=dftb_parameter, cutoff=brcut)
1204 basis_rcut = max(basis_rcut, brcut)
1205 END IF
1206 END IF
1207
1208 IF (PRESENT(total_zeff_corr)) THEN
1209 IF (ASSOCIATED(qs_kind%all_potential)) THEN
1210 CALL get_potential(potential=qs_kind%all_potential, &
1211 zeff=zeff, zeff_correction=zeff_correction)
1212 ELSE IF (ASSOCIATED(qs_kind%gth_potential)) THEN
1213 CALL get_potential(potential=qs_kind%gth_potential, &
1214 zeff=zeff, zeff_correction=zeff_correction)
1215 ELSE IF (ASSOCIATED(qs_kind%sgp_potential)) THEN
1216 CALL get_potential(potential=qs_kind%sgp_potential, &
1217 zeff=zeff, zeff_correction=zeff_correction)
1218 ELSE
1219 zeff = 0.0_dp
1220 zeff_correction = 0.0_dp
1221 END IF
1222 total_zeff_corr = total_zeff_corr + qs_kind_set(ikind)%natom*zeff_correction
1223 END IF
1224
1225 IF (PRESENT(all_potential_present)) THEN
1226 IF (ASSOCIATED(all_potential)) THEN
1227 all_potential_present = .true.
1228 END IF
1229 END IF
1230
1231 IF (PRESENT(tnadd_potential_present)) THEN
1232 IF (ASSOCIATED(tnadd_potential)) THEN
1233 tnadd_potential_present = .true.
1234 END IF
1235 END IF
1236
1237 IF (PRESENT(gth_potential_present)) THEN
1238 IF (ASSOCIATED(gth_potential)) THEN
1239 gth_potential_present = .true.
1240 END IF
1241 END IF
1242
1243 IF (PRESENT(sgp_potential_present)) THEN
1244 IF (ASSOCIATED(sgp_potential)) THEN
1245 sgp_potential_present = .true.
1246 END IF
1247 END IF
1248
1249 IF (PRESENT(cneo_potential_present)) THEN
1250 IF (ASSOCIATED(cneo_potential)) THEN
1251 cneo_potential_present = .true.
1252 END IF
1253 END IF
1254
1255 IF (PRESENT(nkind_q)) THEN
1256 IF (ASSOCIATED(cneo_potential)) THEN
1257 nkind_q = nkind_q + 1
1258 END IF
1259 END IF
1260
1261 IF (PRESENT(natom_q)) THEN
1262 IF (ASSOCIATED(cneo_potential)) THEN
1263 natom_q = natom_q + qs_kind_set(ikind)%natom
1264 END IF
1265 END IF
1266
1267 IF (PRESENT(paw_atom_present)) THEN
1268 IF (paw_atom) THEN
1269 paw_atom_present = .true.
1270 END IF
1271 END IF
1272
1273 IF (PRESENT(dft_plus_u_atom_present)) THEN
1274 IF (dft_plus_u_atom) THEN
1275 dft_plus_u_atom_present = .true.
1276 END IF
1277 END IF
1278
1279 IF (PRESENT(max_ngrid_rad)) THEN
1280 max_ngrid_rad = max(max_ngrid_rad, ngrid_rad)
1281 END IF
1282
1283 IF (PRESENT(max_sph_harm)) THEN
1284 max_sph_harm = max(max_sph_harm, max_s_harm)
1285 END IF
1286
1287 IF (PRESENT(maxg_iso_not0)) THEN
1288 maxg_iso_not0 = max(maxg_iso_not0, max_iso_not0)
1289 END IF
1290
1291 IF (PRESENT(lmax_rho0)) THEN
1292 lmax_rho0 = max(lmax_rho0, lmax_rho0_kind)
1293 END IF
1294
1295 IF (PRESENT(npgf_seg)) THEN
1296 IF (ASSOCIATED(tmp_basis_set)) THEN
1297 CALL get_gto_basis_set(gto_basis_set=tmp_basis_set, npgf_seg_sum=n)
1298 npgf_seg = npgf_seg + n*qs_kind_set(ikind)%natom
1299 END IF
1300 END IF
1301
1302 END DO
1303 ELSE
1304 cpabort("The pointer qs_kind_set is not associated")
1305 END IF
1306
1307 END SUBROUTINE get_qs_kind_set
1308
1309! **************************************************************************************************
1310!> \brief Initialise an atomic kind data set.
1311!> \param qs_kind ...
1312!> \author Creation (11.01.2002,MK)
1313!> 20.09.2002 adapted for pol/kg use, gtb
1314! **************************************************************************************************
1315 SUBROUTINE init_qs_kind(qs_kind)
1316 TYPE(qs_kind_type), POINTER :: qs_kind
1317
1318 CHARACTER(len=*), PARAMETER :: routinen = 'init_qs_kind'
1319
1320 CHARACTER(LEN=default_string_length) :: basis_type
1321 INTEGER :: handle, i
1322 TYPE(gto_basis_set_type), POINTER :: tmp_basis_set
1323
1324 CALL timeset(routinen, handle)
1325
1326 cpassert(ASSOCIATED(qs_kind))
1327
1328 IF (ASSOCIATED(qs_kind%gth_potential)) THEN
1329 CALL init_potential(qs_kind%gth_potential)
1330 ELSEIF (ASSOCIATED(qs_kind%sgp_potential)) THEN
1331 CALL init_potential(qs_kind%sgp_potential)
1332 END IF
1333
1334 DO i = 1, SIZE(qs_kind%basis_sets, 1)
1335 NULLIFY (tmp_basis_set)
1336 CALL get_basis_from_container(qs_kind%basis_sets, basis_set=tmp_basis_set, &
1337 inumbas=i, basis_type=basis_type)
1338 IF (basis_type == "") cycle
1339 IF (basis_type == "AUX") THEN
1340 IF (tmp_basis_set%norm_type < 0) tmp_basis_set%norm_type = 1
1341 CALL init_aux_basis_set(tmp_basis_set)
1342 ELSE
1343 IF (tmp_basis_set%norm_type < 0) tmp_basis_set%norm_type = 2
1344 CALL init_orb_basis_set(tmp_basis_set)
1345 END IF
1346 END DO
1347
1348 CALL timestop(handle)
1349
1350 END SUBROUTINE init_qs_kind
1351
1352! **************************************************************************************************
1353!> \brief Initialise an atomic kind set data set.
1354!> \param qs_kind_set ...
1355!> \author - Creation (17.01.2002,MK)
1356!> - 20.09.2002 para_env passed (gt)
1357! **************************************************************************************************
1358 SUBROUTINE init_qs_kind_set(qs_kind_set)
1359
1360 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
1361
1362 CHARACTER(len=*), PARAMETER :: routinen = 'init_qs_kind_set'
1363
1364 INTEGER :: handle, ikind
1365 TYPE(qs_kind_type), POINTER :: qs_kind
1366
1367 CALL timeset(routinen, handle)
1368
1369 IF (.NOT. ASSOCIATED(qs_kind_set)) THEN
1370 cpabort("init_qs_kind_set: The pointer qs_kind_set is not associated")
1371 END IF
1372
1373 DO ikind = 1, SIZE(qs_kind_set)
1374 qs_kind => qs_kind_set(ikind)
1375 CALL init_qs_kind(qs_kind)
1376 END DO
1377
1378 CALL timestop(handle)
1379
1380 END SUBROUTINE init_qs_kind_set
1381
1382! **************************************************************************************************
1383!> \brief ...
1384!> \param qs_kind_set ...
1385!> \param qs_control ...
1386!> \param force_env_section ...
1387!> \param modify_qs_control whether the qs_control should be modified
1388! **************************************************************************************************
1389 SUBROUTINE init_gapw_basis_set(qs_kind_set, qs_control, force_env_section, modify_qs_control)
1390
1391 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
1392 TYPE(qs_control_type), POINTER :: qs_control
1393 TYPE(section_vals_type), POINTER :: force_env_section
1394 LOGICAL, OPTIONAL :: modify_qs_control
1395
1396 CHARACTER(LEN=default_string_length) :: bsname
1397 INTEGER :: bas1c, ikind, ilevel, nkind
1398 LOGICAL :: gpw, my_mod_control, paw_atom
1399 REAL(dp) :: max_rad_local_type, rc
1400 TYPE(gto_basis_set_type), POINTER :: basis_1c, orb_basis, soft_basis
1401 TYPE(paw_proj_set_type), POINTER :: paw_proj
1402 TYPE(qs_kind_type), POINTER :: qs_kind
1403
1404 my_mod_control = .true.
1405 IF (PRESENT(modify_qs_control)) THEN
1406 my_mod_control = modify_qs_control
1407 END IF
1408
1409 IF (ASSOCIATED(qs_kind_set)) THEN
1410
1411 IF (my_mod_control) qs_control%gapw_control%non_paw_atoms = .false.
1412 nkind = SIZE(qs_kind_set)
1413
1414 DO ikind = 1, nkind
1415
1416 qs_kind => qs_kind_set(ikind)
1417
1418 CALL get_qs_kind(qs_kind=qs_kind, basis_set=orb_basis)
1419 CALL get_qs_kind(qs_kind=qs_kind, hard_radius=rc, &
1420 max_rad_local=max_rad_local_type, gpw_type_forced=gpw)
1421
1422 NULLIFY (soft_basis)
1423 CALL allocate_gto_basis_set(soft_basis)
1424 ! Quantum nuclear wave functions are very localized. Even if soft
1425 ! electronic basis is used, the atomic kind needs PAW treatment
1426 ! because the local Hartree potential is needed.
1427 CALL create_soft_basis(orb_basis, soft_basis, &
1428 qs_control%gapw_control%eps_fit, rc, paw_atom, &
1429 (qs_control%gapw_control%force_paw .OR. &
1430 ASSOCIATED(qs_kind%cneo_potential)), gpw)
1431 CALL add_basis_set_to_container(qs_kind%basis_sets, soft_basis, "ORB_SOFT")
1432 CALL set_qs_kind(qs_kind=qs_kind, paw_atom=paw_atom)
1433
1434 bas1c = qs_control%gapw_control%basis_1c
1435 NULLIFY (basis_1c)
1436 SELECT CASE (bas1c)
1437 CASE (gapw_1c_orb)
1438 ilevel = 0
1439 CASE (gapw_1c_small)
1440 ilevel = 1
1441 CASE (gapw_1c_medium)
1442 ilevel = 2
1443 CASE (gapw_1c_large)
1444 ilevel = 3
1445 CASE (gapw_1c_very_large)
1446 ilevel = 4
1447 CASE DEFAULT
1448 cpabort("basis_1c type")
1449 END SELECT
1450 CALL remove_basis_from_container(qs_kind%basis_sets, basis_type="GAPW_1C")
1451 CALL create_1c_basis(orb_basis, soft_basis, basis_1c, ilevel)
1452 CALL get_gto_basis_set(gto_basis_set=orb_basis, name=bsname)
1453 basis_1c%name = trim(bsname)//"_1c"
1454 CALL add_basis_set_to_container(qs_kind%basis_sets, basis_1c, "GAPW_1C")
1455 IF (paw_atom) THEN
1456 CALL allocate_paw_proj_set(qs_kind%paw_proj_set)
1457 CALL get_qs_kind(qs_kind=qs_kind, paw_proj_set=paw_proj)
1458 CALL projectors(paw_proj, basis_1c, orb_basis, rc, qs_control, &
1459 max_rad_local_type, force_env_section)
1460 ELSE
1461 IF (my_mod_control) qs_control%gapw_control%non_paw_atoms = .true.
1462 END IF
1463
1464 ! grid_atom and harmonics are allocated even if NOT PAW_ATOM
1465 NULLIFY (qs_kind%grid_atom, qs_kind%harmonics)
1466 CALL allocate_grid_atom(qs_kind%grid_atom)
1467 CALL allocate_harmonics_atom(qs_kind%harmonics)
1468
1469 END DO
1470
1471 IF (my_mod_control) THEN
1472 IF (qs_control%gapw_control%non_paw_atoms) THEN
1473 qs_control%gapw_control%nopaw_as_gpw = .true.
1474 ELSE
1475 qs_control%gapw_control%nopaw_as_gpw = .false.
1476 END IF
1477 END IF
1478 ELSE
1479 cpabort("The pointer qs_kind_set is not associated")
1480 END IF
1481
1482 END SUBROUTINE init_gapw_basis_set
1483! **************************************************************************************************
1484!> \brief ...
1485!> \param qs_kind_set ...
1486! **************************************************************************************************
1487 SUBROUTINE init_gapw_nlcc(qs_kind_set)
1488
1489 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
1490
1491 INTEGER :: i, ic, ikind, n_nlcc, nc, nexp_nlcc, &
1492 nkind, nr
1493 INTEGER, DIMENSION(:), POINTER :: nct_nlcc
1494 LOGICAL :: nlcc, nlcc_type, paw_atom
1495 REAL(dp) :: alpha, coa, cval
1496 REAL(kind=dp), DIMENSION(:), POINTER :: a_nlcc, alpha_nlcc, c_nlcc, fe, rc, rr
1497 REAL(kind=dp), DIMENSION(:, :), POINTER :: cval_nlcc, den
1498 TYPE(gth_potential_type), POINTER :: gth_potential
1499 TYPE(qs_kind_type), POINTER :: qs_kind
1500 TYPE(sgp_potential_type), POINTER :: sgp_potential
1501
1502 IF (ASSOCIATED(qs_kind_set)) THEN
1503 nlcc = has_nlcc(qs_kind_set)
1504 IF (nlcc) THEN
1505 nkind = SIZE(qs_kind_set)
1506 DO ikind = 1, nkind
1507 qs_kind => qs_kind_set(ikind)
1508 CALL get_qs_kind(qs_kind, paw_atom=paw_atom)
1509 IF (paw_atom) THEN
1510 CALL get_qs_kind(qs_kind, gth_potential=gth_potential)
1511 CALL get_qs_kind(qs_kind, sgp_potential=sgp_potential)
1512 IF (ASSOCIATED(gth_potential)) THEN
1513 CALL get_potential(potential=gth_potential, nlcc_present=nlcc_type, &
1514 nexp_nlcc=nexp_nlcc, alpha_nlcc=alpha_nlcc, nct_nlcc=nct_nlcc, cval_nlcc=cval_nlcc)
1515 IF (nlcc_type) THEN
1516 nr = qs_kind%grid_atom%nr
1517 rr => qs_kind%grid_atom%rad
1518 ALLOCATE (qs_kind%nlcc_pot(nr, 2), rc(nr), fe(nr))
1519 den => qs_kind%nlcc_pot
1520 den = 0.0_dp
1521 DO i = 1, nexp_nlcc
1522 alpha = alpha_nlcc(i)
1523 rc(:) = rr(:)/alpha
1524 fe(:) = exp(-0.5_dp*rc(:)*rc(:))
1525 nc = nct_nlcc(i)
1526 DO ic = 1, nc
1527 cval = cval_nlcc(ic, i)
1528 coa = cval/alpha
1529 den(:, 1) = den(:, 1) + fe(:)*rc**(2*ic - 2)*cval
1530 den(:, 2) = den(:, 2) - fe(:)*rc**(2*ic - 1)*coa
1531 IF (ic > 1) THEN
1532 den(:, 2) = den(:, 2) + real(2*ic - 2, dp)*fe(:)*rc**(2*ic - 3)*coa
1533 END IF
1534 END DO
1535 END DO
1536 DEALLOCATE (rc, fe)
1537 END IF
1538 ELSE IF (ASSOCIATED(sgp_potential)) THEN
1539 CALL get_potential(potential=sgp_potential, has_nlcc=nlcc_type, &
1540 n_nlcc=n_nlcc, a_nlcc=a_nlcc, c_nlcc=c_nlcc)
1541 IF (nlcc_type) THEN
1542 nr = qs_kind%grid_atom%nr
1543 rr => qs_kind%grid_atom%rad
1544 ALLOCATE (qs_kind%nlcc_pot(nr, 2), rc(nr), fe(nr))
1545 den => qs_kind%nlcc_pot
1546 den = 0.0_dp
1547 DO i = 1, n_nlcc
1548 alpha = a_nlcc(i)
1549 fe(:) = exp(-alpha*rr(:)*rr(:))
1550 cval = c_nlcc(i)
1551 den(:, 1) = den(:, 1) + cval*fe(:)
1552 den(:, 2) = den(:, 2) - 2.0_dp*alpha*cval*rr(:)*fe(:)
1553 END DO
1554 DEALLOCATE (rc, fe)
1555 END IF
1556 ELSE
1557 ! skip
1558 END IF
1559 END IF
1560 END DO
1561 END IF
1562 ELSE
1563 cpabort("The pointer qs_kind_set is not associated")
1564 END IF
1565
1566 END SUBROUTINE init_gapw_nlcc
1567
1568! **************************************************************************************************
1569!> \brief ...
1570!> \param qs_kind_set ...
1571!> \param qs_control ...
1572! **************************************************************************************************
1573 SUBROUTINE init_cneo_basis_set(qs_kind_set, qs_control)
1574
1575 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
1576 TYPE(qs_control_type), POINTER :: qs_control
1577
1578 INTEGER :: ikind, nkind
1579 LOGICAL :: paw_atom
1580 REAL(dp) :: rc
1581 TYPE(gto_basis_set_type), POINTER :: orb_basis, soft_basis
1582 TYPE(qs_kind_type), POINTER :: qs_kind
1583
1584 IF (ASSOCIATED(qs_kind_set)) THEN
1585
1586 nkind = SIZE(qs_kind_set)
1587
1588 DO ikind = 1, nkind
1589
1590 qs_kind => qs_kind_set(ikind)
1591
1592 IF (ASSOCIATED(qs_kind%cneo_potential)) THEN
1593 CALL get_qs_kind(qs_kind=qs_kind, basis_set=orb_basis, basis_type="NUC", &
1594 hard_radius=rc)
1595
1596 NULLIFY (soft_basis)
1597 CALL allocate_gto_basis_set(soft_basis)
1598 CALL create_soft_basis(orb_basis, soft_basis, &
1599 qs_control%gapw_control%eps_fit/ &
1600 sqrt(qs_kind%cneo_potential%zeff), &
1601 rc, paw_atom, .true., .false.)
1602 CALL add_basis_set_to_container(qs_kind%basis_sets, soft_basis, "NUC_SOFT")
1603 END IF
1604
1605 END DO
1606
1607 ELSE
1608 cpabort("The pointer qs_kind_set is not associated")
1609 END IF
1610
1611 END SUBROUTINE init_cneo_basis_set
1612
1613! **************************************************************************************************
1614!> \brief Read an atomic kind data set from the input file.
1615!> \param qs_kind ...
1616!> \param kind_section ...
1617!> \param para_env ...
1618!> \param force_env_section ...
1619!> \param no_fail ...
1620!> \param method_id ...
1621!> \param silent ...
1622!> \par History
1623!> - Creation (09.02.2002,MK)
1624!> - 20.09.2002,gt: adapted for POL/KG use (elp_potential)
1625!> - 05.03.2010: split elp_potential into fist_potential and kg_potential
1626! **************************************************************************************************
1627 SUBROUTINE read_qs_kind(qs_kind, kind_section, para_env, force_env_section, &
1628 no_fail, method_id, silent)
1629
1630 TYPE(qs_kind_type), INTENT(INOUT) :: qs_kind
1631 TYPE(section_vals_type), POINTER :: kind_section
1632 TYPE(mp_para_env_type), POINTER :: para_env
1633 TYPE(section_vals_type), POINTER :: force_env_section
1634 LOGICAL, INTENT(IN) :: no_fail
1635 INTEGER, INTENT(IN) :: method_id
1636 LOGICAL, INTENT(IN) :: silent
1637
1638 CHARACTER(LEN=*), PARAMETER :: routinen = 'read_qs_kind'
1639 INTEGER, PARAMETER :: maxbas = 20
1640
1641 CHARACTER(LEN=2) :: element_symbol
1642 CHARACTER(len=default_path_length) :: kg_potential_fn_kind, &
1643 potential_file_name, potential_fn_kind
1644 CHARACTER(LEN=default_string_length) :: akind_name, basis_type, keyword, &
1645 kgpot_name, kgpot_type, &
1646 potential_name, potential_type, tmp
1647 CHARACTER(LEN=default_string_length), DIMENSION(4) :: description
1648 CHARACTER(LEN=default_string_length), &
1649 DIMENSION(:), POINTER :: tmpstringlist
1650 CHARACTER(LEN=default_string_length), &
1651 DIMENSION(maxbas) :: basis_set_form, basis_set_name, &
1652 basis_set_type
1653 INTEGER :: handle, i, i_rep, iounit, ipaodesc, ipaopot, ipos, j, jj, k_rep, l, m, n_rep, &
1654 nb_rep, nexp, ngauss, nlcc, nloc, nnl, norbitals, npaodesc, npaopot, nppnl, nspin, nu, z
1655 INTEGER, DIMENSION(:), POINTER :: add_el, elec_conf, orbitals
1656 LOGICAL :: check, ecp_semi_local, explicit, explicit_basis, explicit_j, explicit_kgpot, &
1657 explicit_potential, explicit_u, explicit_u_m_j, nobasis, nobasis_nuc, section_enabled, &
1658 subsection_enabled, update_input
1659 REAL(kind=dp) :: alpha, ccore, mass, r, rc, &
1660 zeff_correction
1661 REAL(kind=dp), DIMENSION(6) :: error
1662 REAL(kind=dp), DIMENSION(:), POINTER :: a_nl, aloc, anlcc, cloc, cnlcc, nelec
1663 REAL(kind=dp), DIMENSION(:, :), POINTER :: h_nl
1664 REAL(kind=dp), DIMENSION(:, :, :), POINTER :: c_nl
1665 TYPE(atom_ecppot_type) :: ecppot
1666 TYPE(atom_sgp_potential_type) :: sgppot
1667 TYPE(atom_upfpot_type) :: upfpot
1668 TYPE(cp_logger_type), POINTER :: logger
1669 TYPE(gto_basis_set_type), POINTER :: orb_basis_set, sup_basis_set, &
1670 tmp_basis_set
1671 TYPE(section_vals_type), POINTER :: basis_section, bs_section, dft_plus_u_section, &
1672 dft_section, enforce_occupation_section, kgpot_section, pao_desc_section, &
1673 pao_pot_section, potential_section, spin_section
1674 TYPE(sto_basis_set_type), POINTER :: sto_basis_set
1675
1676 CALL timeset(routinen, handle)
1677
1678 NULLIFY (logger)
1679 logger => cp_get_default_logger()
1680 iounit = cp_logger_get_default_io_unit(logger)
1681
1682 NULLIFY (elec_conf)
1683
1684 update_input = .true.
1685 basis_set_name(:) = ""
1686 basis_set_type(:) = ""
1687 basis_set_form(:) = ""
1688 potential_name = ""
1689 potential_type = ""
1690 kgpot_name = ""
1691 kgpot_type = ""
1692 z = -1
1693 zeff_correction = 0.0_dp
1694 explicit = .false.
1695 explicit_basis = .false.
1696 explicit_j = .false.
1697 explicit_kgpot = .false.
1698 explicit_potential = .false.
1699 explicit_u = .false.
1700 explicit_u_m_j = .false.
1701
1702 dft_section => section_vals_get_subs_vals(force_env_section, "DFT")
1703 CALL section_vals_get(kind_section, n_repetition=n_rep)
1704 k_rep = -1
1705 akind_name = qs_kind%name
1706 CALL uppercase(akind_name)
1707 ! First we use the atom_name to find out the proper KIND section
1708 DO i_rep = 1, n_rep
1709 CALL section_vals_val_get(kind_section, "_SECTION_PARAMETERS_", &
1710 c_val=keyword, i_rep_section=i_rep)
1711 CALL uppercase(keyword)
1712 IF (keyword == akind_name) THEN
1713 k_rep = i_rep
1714 EXIT
1715 END IF
1716 END DO
1717 ! The search for the KIND section failed.. check for a QM/MM link atom
1718 IF (k_rep < 1) THEN
1719 ipos = index(qs_kind%name, "_")
1720 IF (((ipos == 2) .OR. (ipos == 3)) .AND. (index(qs_kind%name, "_ghost") == 0)) THEN
1721 ! If the atm_name could not match any KIND section it maybe be a QM/MM link atom.
1722 ! ghost atoms will be treated differently.
1723 akind_name = qs_kind%name(1:ipos - 1)
1724 CALL uppercase(akind_name)
1725 DO i_rep = 1, n_rep
1726 CALL section_vals_val_get(kind_section, "_SECTION_PARAMETERS_", &
1727 c_val=keyword, i_rep_section=i_rep)
1728 CALL uppercase(keyword)
1729 IF (keyword == akind_name) THEN
1730 k_rep = i_rep
1731 EXIT
1732 END IF
1733 END DO
1734 END IF
1735 END IF
1736 ! The search for the KIND section failed.. check element_symbol
1737 IF (k_rep < 1) THEN
1738 ! If it's not a link atom let's check for the element and map
1739 ! the KIND section to the element.
1740 element_symbol = qs_kind%element_symbol(1:2)
1741 CALL uppercase(element_symbol)
1742 DO i_rep = 1, n_rep
1743 CALL section_vals_val_get(kind_section, "_SECTION_PARAMETERS_", &
1744 c_val=keyword, i_rep_section=i_rep)
1745 CALL uppercase(keyword)
1746 IF (keyword == element_symbol) THEN
1747 k_rep = i_rep
1748 EXIT
1749 END IF
1750 END DO
1751 END IF
1752 ! In case it should not really match any possible KIND section
1753 ! let's look if a default one is defined..
1754 IF (k_rep < 1) THEN
1755 DO i_rep = 1, n_rep
1756 CALL section_vals_val_get(kind_section, "_SECTION_PARAMETERS_", &
1757 c_val=keyword, i_rep_section=i_rep)
1758 CALL uppercase(keyword)
1759 IF (keyword == "DEFAULT") THEN
1760 update_input = .false.
1761 k_rep = i_rep
1762 EXIT
1763 END IF
1764 END DO
1765 END IF
1766 IF (k_rep < 0 .AND. (.NOT. no_fail)) THEN
1767 CALL cp_abort(__location__, &
1768 "No &KIND section was possible to associate to the atomic kind <"// &
1769 trim(akind_name)//">. The KIND section were also scanned for the"// &
1770 " corresponding element <"//trim(qs_kind%element_symbol)//">"// &
1771 " and for the DEFAULT section but no match was found. Check your input file!")
1772 END IF
1773 ! Retrieve information on element
1774 CALL get_ptable_info(qs_kind%element_symbol, ielement=z)
1775
1776 ! Normal parsing of the KIND section
1777 IF (k_rep > 0) THEN
1778 ! new style basis set input
1779 CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
1780 keyword_name="BASIS_SET", &
1781 explicit=explicit, &
1782 n_rep_val=nb_rep)
1783 IF (.NOT. explicit) nb_rep = 0
1784 cpassert(nb_rep <= maxbas)
1785 DO i = 1, nb_rep
1786 CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
1787 keyword_name="BASIS_SET", i_rep_val=i, c_vals=tmpstringlist)
1788 IF (SIZE(tmpstringlist) == 1) THEN
1789 ! default is orbital type and GTO
1790 basis_set_type(i) = "ORB"
1791 basis_set_form(i) = "GTO"
1792 basis_set_name(i) = tmpstringlist(1)
1793 ELSEIF (SIZE(tmpstringlist) == 2) THEN
1794 ! default is GTO
1795 basis_set_type(i) = tmpstringlist(1)
1796 basis_set_form(i) = "GTO"
1797 basis_set_name(i) = tmpstringlist(2)
1798 ELSEIF (SIZE(tmpstringlist) == 3) THEN
1799 basis_set_type(i) = tmpstringlist(1)
1800 basis_set_form(i) = tmpstringlist(2)
1801 basis_set_name(i) = tmpstringlist(3)
1802 ELSE
1803 CALL cp_abort(__location__, &
1804 "invalid number of BASIS_SET keyword parameters: BASIS_SET [<TYPE>] [<FORM>] <NAME>")
1805 END IF
1806 ! check that we have a valid basis set form
1807 IF (basis_set_form(i) /= "GTO" .AND. basis_set_form(i) /= "STO") THEN
1808 cpabort("invalid BASIS_SET FORM parameter")
1809 END IF
1810 END DO
1811
1812 ! parse PAO keywords
1813 CALL section_vals_val_get(kind_section, keyword_name="PAO_BASIS_SIZE", i_rep_section=k_rep, &
1814 i_val=qs_kind%pao_basis_size)
1815 CALL section_vals_val_get(kind_section, keyword_name="PAO_MODEL_FILE", i_rep_section=k_rep, &
1816 explicit=explicit)
1817 IF (explicit) THEN
1818 CALL section_vals_val_get(kind_section, keyword_name="PAO_MODEL_FILE", i_rep_section=k_rep, &
1819 c_val=qs_kind%pao_model_file)
1820 END IF
1821
1822 ! parse PAO_POTENTIAL sections
1823 pao_pot_section => section_vals_get_subs_vals(kind_section, "PAO_POTENTIAL", i_rep_section=k_rep)
1824 CALL section_vals_get(pao_pot_section, n_repetition=npaopot)
1825 ALLOCATE (qs_kind%pao_potentials(npaopot))
1826 DO ipaopot = 1, npaopot
1827 CALL section_vals_val_get(pao_pot_section, keyword_name="MAXL", i_rep_section=ipaopot, &
1828 i_val=qs_kind%pao_potentials(ipaopot)%maxl)
1829 CALL section_vals_val_get(pao_pot_section, keyword_name="MAX_PROJECTOR", i_rep_section=ipaopot, &
1830 i_val=qs_kind%pao_potentials(ipaopot)%max_projector)
1831 CALL section_vals_val_get(pao_pot_section, keyword_name="BETA", i_rep_section=ipaopot, &
1832 r_val=qs_kind%pao_potentials(ipaopot)%beta)
1833 CALL section_vals_val_get(pao_pot_section, keyword_name="WEIGHT", i_rep_section=ipaopot, &
1834 r_val=qs_kind%pao_potentials(ipaopot)%weight)
1835 END DO
1836
1837 ! parse PAO_DESCRIPTOR sections
1838 pao_desc_section => section_vals_get_subs_vals(kind_section, "PAO_DESCRIPTOR", i_rep_section=k_rep)
1839 CALL section_vals_get(pao_desc_section, n_repetition=npaodesc)
1840 ALLOCATE (qs_kind%pao_descriptors(npaodesc))
1841 DO ipaodesc = 1, npaodesc
1842 CALL section_vals_val_get(pao_desc_section, keyword_name="BETA", i_rep_section=ipaodesc, &
1843 r_val=qs_kind%pao_descriptors(ipaodesc)%beta)
1844 CALL section_vals_val_get(pao_desc_section, keyword_name="SCREENING", i_rep_section=ipaodesc, &
1845 r_val=qs_kind%pao_descriptors(ipaodesc)%screening)
1846 CALL section_vals_val_get(pao_desc_section, keyword_name="WEIGHT", i_rep_section=ipaodesc, &
1847 r_val=qs_kind%pao_descriptors(ipaodesc)%weight)
1848 END DO
1849
1850 ! parse ELEC_CONF
1851 CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
1852 keyword_name="ELEC_CONF", n_rep_val=i)
1853 IF (i > 0) THEN
1854 CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
1855 keyword_name="ELEC_CONF", i_vals=elec_conf)
1856 CALL set_qs_kind(qs_kind, elec_conf=elec_conf)
1857 END IF
1858 CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
1859 keyword_name="CORE_CORRECTION", r_val=zeff_correction)
1860 ! parse POTENTIAL
1861 CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
1862 keyword_name="POTENTIAL_FILE_NAME", c_val=potential_fn_kind)
1863 CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
1864 keyword_name="POTENTIAL_TYPE", c_val=potential_type)
1865 CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
1866 explicit=explicit, keyword_name="POTENTIAL", c_vals=tmpstringlist)
1867 IF (explicit) THEN
1868 IF (SIZE(tmpstringlist) == 1) THEN
1869 ! old type of input: start of name defines type
1870 potential_name = tmpstringlist(1)
1871 IF (potential_type == "") THEN
1872 ipos = index(potential_name, "-")
1873 IF (ipos > 1) THEN
1874 potential_type = potential_name(:ipos - 1)
1875 ELSE
1876 potential_type = potential_name
1877 END IF
1878 END IF
1879 ELSEIF (SIZE(tmpstringlist) == 2) THEN
1880 potential_type = tmpstringlist(1)
1881 potential_name = tmpstringlist(2)
1882 ELSE
1883 cpabort("POTENTIAL input list is not correct")
1884 END IF
1885 END IF
1886 CALL uppercase(potential_type)
1887
1888 ! Parse KG POTENTIAL
1889 CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
1890 keyword_name="KG_POTENTIAL_FILE_NAME", c_val=kg_potential_fn_kind)
1891 CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
1892 keyword_name="KG_POTENTIAL", c_val=kgpot_name)
1893
1894 ! Semi-local vs. full nonlocal form of ECPs
1895 CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
1896 keyword_name="ECP_SEMI_LOCAL", l_val=ecp_semi_local)
1897
1898 ! Assign atomic covalent radius
1899 qs_kind%covalent_radius = ptable(z)%covalent_radius*bohr
1900 CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
1901 keyword_name="COVALENT_RADIUS", r_val=r)
1902 IF (r > 0.0_dp) qs_kind%covalent_radius = r
1903
1904 ! Assign atomic van der Waals radius
1905 qs_kind%vdw_radius = ptable(z)%vdw_radius*bohr
1906 CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
1907 keyword_name="VDW_RADIUS", r_val=r)
1908 IF (r > 0.0_dp) qs_kind%vdw_radius = r
1909
1910 ! Assign atom dependent defaults, only H special case
1911 CALL section_vals_val_get(kind_section, i_rep_section=k_rep, n_rep_val=i, &
1912 keyword_name="HARD_EXP_RADIUS")
1913 IF (i == 0) THEN
1914 IF (z == 1) THEN
1915 qs_kind%hard_radius = 1.2_dp
1916 ELSE
1917 qs_kind%hard_radius = 0.8_dp*bohr
1918 END IF
1919 ELSE
1920 CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
1921 keyword_name="HARD_EXP_RADIUS", r_val=qs_kind%hard_radius)
1922 END IF
1923
1924 ! assign atom dependent defaults, only H special case
1925 CALL section_vals_val_get(kind_section, i_rep_section=k_rep, n_rep_val=i, &
1926 keyword_name="RHO0_EXP_RADIUS")
1927 IF (i == 0) THEN
1928 qs_kind%hard0_radius = qs_kind%hard_radius
1929 ELSE
1930 CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
1931 keyword_name="RHO0_EXP_RADIUS", r_val=qs_kind%hard0_radius)
1932 END IF
1933 IF (qs_kind%hard_radius < qs_kind%hard0_radius) &
1934 cpabort("rc0 should be <= rc")
1935
1936 CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
1937 keyword_name="MAX_RAD_LOCAL", r_val=qs_kind%max_rad_local)
1938 CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
1939 keyword_name="LEBEDEV_GRID", i_val=qs_kind%ngrid_ang)
1940 IF (qs_kind%ngrid_ang <= 0) &
1941 cpabort("# point lebedev grid < 0")
1942 CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
1943 keyword_name="RADIAL_GRID", i_val=qs_kind%ngrid_rad)
1944 IF (qs_kind%ngrid_rad <= 0) &
1945 cpabort("# point radial grid < 0")
1946 CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
1947 keyword_name="GPW_TYPE", l_val=qs_kind%gpw_type_forced)
1948 CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
1949 keyword_name="GHOST", l_val=qs_kind%ghost)
1950 CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
1951 keyword_name="FLOATING_BASIS_CENTER", l_val=qs_kind%floating)
1952 CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
1953 keyword_name="NO_OPTIMIZE", l_val=qs_kind%no_optimize)
1954
1955 ! Magnetization
1956 CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
1957 keyword_name="MAGNETIZATION", r_val=qs_kind%magnetization)
1958 ! DFTB3 param
1959 CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
1960 keyword_name="DFTB3_PARAM", r_val=qs_kind%dudq_dftb3)
1961 CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
1962 keyword_name="LMAX_DFTB", i_val=qs_kind%lmax_dftb)
1963
1964 ! MAOS
1965 CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
1966 keyword_name="MAO", i_val=qs_kind%mao)
1967
1968 ! Use monovalent pseudopotential
1969 CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
1970 keyword_name="MONOVALENT", l_val=qs_kind%monovalent)
1971 IF (qs_kind%monovalent .AND. trim(potential_type) /= 'GTH') THEN
1972 cpabort("Monovalent pseudopotentials currently implemented for GTH only!")
1973 END IF
1974
1975 ! Read the BS subsection of the current atomic kind, if enabled
1976 NULLIFY (bs_section)
1977 bs_section => section_vals_get_subs_vals(kind_section, "BS", &
1978 i_rep_section=k_rep)
1979 section_enabled = .false.
1980 CALL section_vals_val_get(bs_section, "_SECTION_PARAMETERS_", &
1981 l_val=section_enabled)
1982 IF (section_enabled) THEN
1983 ! test for conflict with magnetization
1984 IF (qs_kind%magnetization /= 0.0_dp) THEN
1985 CALL cp_abort(__location__, "BS Section is in conflict with non-zero magnetization "// &
1986 "for this atom kind.")
1987 END IF
1988 qs_kind%bs_occupation = .true.
1989 !Alpha spin
1990 NULLIFY (spin_section)
1991 spin_section => section_vals_get_subs_vals(bs_section, "ALPHA")
1992 CALL section_vals_get(spin_section, explicit=explicit)
1993 IF (explicit) THEN
1994 NULLIFY (add_el)
1995 CALL section_vals_val_get(spin_section, &
1996 keyword_name="NEL", i_vals=add_el)
1997 cpassert(ASSOCIATED(add_el))
1998 ALLOCATE (qs_kind%addel(SIZE(add_el), 2))
1999 qs_kind%addel = 0
2000 qs_kind%addel(1:SIZE(add_el), 1) = add_el(1:SIZE(add_el))
2001 NULLIFY (add_el)
2002 CALL section_vals_val_get(spin_section, &
2003 keyword_name="L", i_vals=add_el)
2004 cpassert(ASSOCIATED(add_el))
2005 cpassert(SIZE(add_el) == SIZE(qs_kind%addel, 1))
2006 ALLOCATE (qs_kind%laddel(SIZE(add_el), 2))
2007 qs_kind%laddel = 0
2008 qs_kind%laddel(1:SIZE(add_el), 1) = add_el(1:SIZE(add_el))
2009 ALLOCATE (qs_kind%naddel(SIZE(add_el), 2))
2010 qs_kind%naddel = 0
2011 NULLIFY (add_el)
2012 CALL section_vals_val_get(spin_section, &
2013 keyword_name="N", n_rep_val=i)
2014 IF (i > 0) THEN
2015 CALL section_vals_val_get(spin_section, &
2016 keyword_name="N", i_vals=add_el)
2017 IF (SIZE(add_el) == SIZE(qs_kind%addel, 1)) THEN
2018 qs_kind%naddel(1:SIZE(add_el), 1) = add_el(1:SIZE(add_el))
2019 END IF
2020 END IF
2021 END IF
2022 ! Beta spin
2023 NULLIFY (spin_section)
2024 spin_section => section_vals_get_subs_vals(bs_section, "BETA")
2025 CALL section_vals_get(spin_section, explicit=explicit)
2026 IF (explicit) THEN
2027 NULLIFY (add_el)
2028 CALL section_vals_val_get(spin_section, &
2029 keyword_name="NEL", i_vals=add_el)
2030 cpassert(SIZE(add_el) == SIZE(qs_kind%addel, 1))
2031 qs_kind%addel(1:SIZE(add_el), 2) = add_el(1:SIZE(add_el))
2032 qs_kind%addel(:, :) = qs_kind%addel(:, :)
2033 NULLIFY (add_el)
2034 CALL section_vals_val_get(spin_section, &
2035 keyword_name="L", i_vals=add_el)
2036 cpassert(SIZE(add_el) == SIZE(qs_kind%addel, 1))
2037 qs_kind%laddel(1:SIZE(add_el), 2) = add_el(1:SIZE(add_el))
2038
2039 CALL section_vals_val_get(spin_section, &
2040 keyword_name="N", n_rep_val=i)
2041 IF (i > 0) THEN
2042 NULLIFY (add_el)
2043 CALL section_vals_val_get(spin_section, &
2044 keyword_name="N", i_vals=add_el)
2045 IF (SIZE(add_el) == SIZE(qs_kind%addel, 1)) THEN
2046 qs_kind%naddel(1:SIZE(add_el), 2) = add_el(1:SIZE(add_el))
2047 END IF
2048 END IF
2049 END IF
2050 END IF
2051
2052 ! Read the DFT+U subsection of the current atomic kind, if enabled
2053
2054 NULLIFY (dft_plus_u_section)
2055 dft_plus_u_section => section_vals_get_subs_vals(kind_section, &
2056 subsection_name="DFT_PLUS_U", &
2057 i_rep_section=k_rep)
2058 section_enabled = .false.
2059 CALL section_vals_val_get(dft_plus_u_section, &
2060 keyword_name="_SECTION_PARAMETERS_", &
2061 l_val=section_enabled)
2062 IF (section_enabled) THEN
2063 ALLOCATE (qs_kind%dft_plus_u)
2064 NULLIFY (qs_kind%dft_plus_u%nelec)
2065 NULLIFY (qs_kind%dft_plus_u%orbitals)
2066 CALL section_vals_val_get(dft_plus_u_section, &
2067 keyword_name="L", &
2068 i_val=l)
2069 qs_kind%dft_plus_u%l = l
2070#if defined(__SIRIUS)
2071 CALL section_vals_val_get(dft_plus_u_section, &
2072 keyword_name="N", &
2073 i_val=nu)
2074 qs_kind%dft_plus_u%n = nu
2075
2076 CALL section_vals_val_get(dft_plus_u_section, &
2077 keyword_name="U", &
2078 r_val=qs_kind%dft_plus_u%U, &
2079 explicit=explicit_u)
2080
2081 CALL section_vals_val_get(dft_plus_u_section, &
2082 keyword_name="J", &
2083 r_val=qs_kind%dft_plus_u%J, &
2084 explicit=explicit_j)
2085
2086 CALL section_vals_val_get(dft_plus_u_section, &
2087 keyword_name="alpha", &
2088 r_val=qs_kind%dft_plus_u%alpha)
2089
2090 CALL section_vals_val_get(dft_plus_u_section, &
2091 keyword_name="beta", &
2092 r_val=qs_kind%dft_plus_u%beta)
2093
2094 CALL section_vals_val_get(dft_plus_u_section, &
2095 keyword_name="J0", &
2096 r_val=qs_kind%dft_plus_u%J0)
2097
2098 CALL section_vals_val_get(dft_plus_u_section, &
2099 keyword_name="occupation", &
2100 r_val=qs_kind%dft_plus_u%occupation)
2101#else
2102 nu = 0
2103#endif
2104
2105 CALL section_vals_val_get(dft_plus_u_section, &
2106 keyword_name="U_MINUS_J", &
2107 r_val=qs_kind%dft_plus_u%u_minus_j_target, &
2108 explicit=explicit_u_m_j)
2109
2110 IF ((explicit_u .OR. explicit_j) .AND. explicit_u_m_j) THEN
2111 cpabort("DFT+U| specifying U or J and U_MINUS_J parameters are mutually exclusive.")
2112 END IF
2113
2114 CALL section_vals_val_get(dft_plus_u_section, &
2115 keyword_name="U_RAMPING", &
2116 r_val=qs_kind%dft_plus_u%u_ramping)
2117 CALL section_vals_val_get(dft_plus_u_section, &
2118 keyword_name="INIT_U_RAMPING_EACH_SCF", &
2119 l_val=qs_kind%dft_plus_u%init_u_ramping_each_scf)
2120 IF (qs_kind%dft_plus_u%u_ramping > 0.0_dp) THEN
2121 qs_kind%dft_plus_u%u_minus_j = 0.0_dp
2122 ELSE
2123 qs_kind%dft_plus_u%u_minus_j = qs_kind%dft_plus_u%u_minus_j_target
2124 END IF
2125 CALL section_vals_val_get(dft_plus_u_section, &
2126 keyword_name="EPS_U_RAMPING", &
2127 r_val=qs_kind%dft_plus_u%eps_u_ramping)
2128
2129 NULLIFY (enforce_occupation_section)
2130 enforce_occupation_section => section_vals_get_subs_vals(dft_plus_u_section, &
2131 subsection_name="ENFORCE_OCCUPATION")
2132 subsection_enabled = .false.
2133 CALL section_vals_val_get(enforce_occupation_section, &
2134 keyword_name="_SECTION_PARAMETERS_", &
2135 l_val=subsection_enabled)
2136 IF (subsection_enabled) THEN
2137 NULLIFY (nelec)
2138 CALL section_vals_val_get(enforce_occupation_section, &
2139 keyword_name="NELEC", &
2140 r_vals=nelec)
2141 nspin = SIZE(nelec)
2142 ALLOCATE (qs_kind%dft_plus_u%nelec(nspin))
2143 qs_kind%dft_plus_u%nelec(:) = nelec(:)
2144 NULLIFY (orbitals)
2145 CALL section_vals_val_get(enforce_occupation_section, &
2146 keyword_name="ORBITALS", &
2147 i_vals=orbitals)
2148 norbitals = SIZE(orbitals)
2149 IF (norbitals <= 0 .OR. norbitals > 2*l + 1) &
2150 CALL cp_abort(__location__, "DFT+U| Invalid number of ORBITALS specified: "// &
2151 "1 to 2*L+1 integer numbers are expected")
2152 ALLOCATE (qs_kind%dft_plus_u%orbitals(norbitals))
2153 qs_kind%dft_plus_u%orbitals(:) = orbitals(:)
2154 NULLIFY (orbitals)
2155 DO m = 1, norbitals
2156 IF (qs_kind%dft_plus_u%orbitals(m) > l) &
2157 cpabort("DFT+U| Invalid orbital magnetic quantum number specified: m > l")
2158 IF (qs_kind%dft_plus_u%orbitals(m) < -l) &
2159 cpabort("DFT+U| Invalid orbital magnetic quantum number specified: m < -l")
2160 DO j = 1, norbitals
2161 IF (j /= m) THEN
2162 IF (qs_kind%dft_plus_u%orbitals(j) == qs_kind%dft_plus_u%orbitals(m)) &
2163 cpabort("DFT+U| An orbital magnetic quantum number was specified twice")
2164 END IF
2165 END DO
2166 END DO
2167 CALL section_vals_val_get(enforce_occupation_section, &
2168 keyword_name="EPS_SCF", &
2169 r_val=qs_kind%dft_plus_u%eps_scf)
2170 CALL section_vals_val_get(enforce_occupation_section, &
2171 keyword_name="MAX_SCF", &
2172 i_val=i)
2173 qs_kind%dft_plus_u%max_scf = max(-1, i)
2174 CALL section_vals_val_get(enforce_occupation_section, &
2175 keyword_name="SMEAR", &
2176 l_val=qs_kind%dft_plus_u%smear)
2177 END IF ! subsection enabled
2178 END IF ! section enabled
2179
2180 END IF
2181
2182 ! Allocate and initialise the orbital basis set data set structure
2183 CALL init_orbital_pointers(5) ! debug the SUN optimizer
2184
2185 ! BASIS and POTENTIAL read only when strictly necessary otherwise, even if not used
2186 ! we just print misleading informations
2187 explicit_basis = .false.
2188 IF (k_rep > 0) THEN
2189 basis_section => section_vals_get_subs_vals(kind_section, "BASIS", i_rep_section=k_rep, &
2190 can_return_null=.true.)
2191 CALL section_vals_get(basis_section, explicit=explicit_basis)
2192 END IF
2193
2194 explicit_potential = .false.
2195 IF (k_rep > 0) THEN
2196 potential_section => section_vals_get_subs_vals(kind_section, "POTENTIAL", &
2197 i_rep_section=k_rep, can_return_null=.true.)
2198 CALL section_vals_get(potential_section, explicit=explicit_potential)
2199 END IF
2200
2201 explicit_kgpot = .false.
2202 IF (k_rep > 0) THEN
2203 kgpot_section => section_vals_get_subs_vals(kind_section, "KG_POTENTIAL", &
2204 i_rep_section=k_rep, can_return_null=.true.)
2205 CALL section_vals_get(kgpot_section, explicit=explicit_kgpot)
2206 END IF
2207
2208 SELECT CASE (method_id)
2211 ! Allocate all_potential
2212 CALL allocate_potential(qs_kind%all_potential)
2213 CALL set_default_all_potential(qs_kind%all_potential, z, zeff_correction)
2214 CALL get_qs_kind(qs_kind, elec_conf=elec_conf)
2215 IF (.NOT. ASSOCIATED(elec_conf)) THEN
2216 CALL get_potential(potential=qs_kind%all_potential, elec_conf=elec_conf)
2217 CALL set_qs_kind(qs_kind, elec_conf=elec_conf)
2218 END IF
2219 cpassert(.NOT. qs_kind%floating)
2220 IF (qs_kind%ghost) THEN
2221 CALL get_qs_kind(qs_kind=qs_kind, elec_conf=elec_conf)
2222 elec_conf(:) = 0
2223 CALL get_potential(potential=qs_kind%all_potential, &
2224 elec_conf=elec_conf)
2225 elec_conf(:) = 0
2226 CALL set_potential(potential=qs_kind%all_potential, &
2227 zeff=0.0_dp, &
2228 zeff_correction=0.0_dp)
2229 END IF
2230
2231 ! Basis set (Parameters)
2232 ! Setup proper semiempirical parameters
2233 check = .NOT. ASSOCIATED(qs_kind%se_parameter)
2234 cpassert(check)
2235 CALL semi_empirical_create(qs_kind%se_parameter)
2236 ! Check if we allow p-orbitals on H
2237 SELECT CASE (z)
2238 CASE (1)
2239 IF (k_rep > 0) THEN
2240 CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
2241 keyword_name="SE_P_ORBITALS_ON_H", l_val=qs_kind%se_parameter%p_orbitals_on_h)
2242 END IF
2243 CASE DEFAULT
2244 ! No special cases for other elements..
2245 END SELECT
2246 ! Set default parameters
2247 CALL section_vals_val_get(dft_section, "QS%SE%STO_NG", i_val=ngauss)
2248 CALL se_param_set_default(qs_kind%se_parameter, z, method_id)
2249 NULLIFY (tmp_basis_set)
2250 CALL init_se_param(qs_kind%se_parameter, tmp_basis_set, ngauss)
2251 CALL add_basis_set_to_container(qs_kind%basis_sets, tmp_basis_set, "ORB")
2252 CALL init_potential(qs_kind%all_potential, itype="BARE", &
2253 zeff=qs_kind%se_parameter%zeff, zeff_correction=zeff_correction)
2254 qs_kind%se_parameter%zeff = qs_kind%se_parameter%zeff - zeff_correction
2255
2256 check = ((potential_name /= '') .OR. explicit_potential) .AND. .NOT. silent
2257 IF (check) &
2258 CALL cp_warn(__location__, &
2259 "Information provided in the input file regarding POTENTIAL for KIND <"// &
2260 trim(qs_kind%name)//"> will be ignored!")
2261
2262 check = ((k_rep > 0) .OR. explicit_basis) .AND. .NOT. silent
2263 IF (check) &
2264 CALL cp_warn(__location__, &
2265 "Information provided in the input file regarding BASIS for KIND <"// &
2266 trim(qs_kind%name)//"> will be ignored!")
2267
2268 CASE (do_method_dftb)
2269 ! Allocate all_potential
2270 CALL allocate_potential(qs_kind%all_potential)
2271 CALL set_default_all_potential(qs_kind%all_potential, z, zeff_correction)
2272 CALL get_qs_kind(qs_kind, elec_conf=elec_conf)
2273 IF (.NOT. ASSOCIATED(elec_conf)) THEN
2274 CALL get_potential(potential=qs_kind%all_potential, elec_conf=elec_conf)
2275 CALL set_qs_kind(qs_kind, elec_conf=elec_conf)
2276 END IF
2277 cpassert(.NOT. qs_kind%floating)
2278 IF (qs_kind%ghost) THEN
2279 CALL get_qs_kind(qs_kind=qs_kind, elec_conf=elec_conf)
2280 elec_conf(:) = 0
2281 CALL get_potential(potential=qs_kind%all_potential, &
2282 elec_conf=elec_conf)
2283 elec_conf(:) = 0
2284 CALL set_potential(potential=qs_kind%all_potential, &
2285 zeff=0.0_dp, &
2286 zeff_correction=0.0_dp)
2287 END IF
2288
2289 check = ((potential_name /= '') .OR. explicit_potential) .AND. .NOT. silent
2290 IF (check) &
2291 CALL cp_warn(__location__, &
2292 "Information provided in the input file regarding POTENTIAL for KIND <"// &
2293 trim(qs_kind%name)//"> will be ignored!")
2294
2295 check = ((k_rep > 0) .OR. explicit_basis) .AND. .NOT. silent
2296 IF (check) &
2297 CALL cp_warn(__location__, &
2298 "Information provided in the input file regarding BASIS for KIND <"// &
2299 trim(qs_kind%name)//"> will be ignored!")
2300
2301 CASE (do_method_xtb)
2302 ! Allocate all_potential
2303 CALL allocate_potential(qs_kind%all_potential)
2304 CALL set_default_all_potential(qs_kind%all_potential, z, zeff_correction)
2305 CALL get_qs_kind(qs_kind, elec_conf=elec_conf)
2306 IF (.NOT. ASSOCIATED(elec_conf)) THEN
2307 CALL get_potential(potential=qs_kind%all_potential, elec_conf=elec_conf)
2308 CALL set_qs_kind(qs_kind, elec_conf=elec_conf)
2309 END IF
2310 cpassert(.NOT. qs_kind%floating)
2311 IF (qs_kind%ghost) THEN
2312 CALL get_qs_kind(qs_kind=qs_kind, elec_conf=elec_conf)
2313 elec_conf(:) = 0
2314 CALL get_potential(potential=qs_kind%all_potential, &
2315 elec_conf=elec_conf)
2316 elec_conf(:) = 0
2317 CALL set_potential(potential=qs_kind%all_potential, &
2318 zeff=0.0_dp, &
2319 zeff_correction=0.0_dp)
2320 END IF
2321
2322 check = ((potential_name /= '') .OR. explicit_potential) .AND. .NOT. silent
2323 IF (check) &
2324 CALL cp_warn(__location__, &
2325 "Information provided in the input file regarding POTENTIAL for KIND <"// &
2326 trim(qs_kind%name)//"> will be ignored!")
2327
2328 check = ((k_rep > 0) .OR. explicit_basis) .AND. .NOT. silent
2329 IF (check) &
2330 CALL cp_warn(__location__, &
2331 "Information provided in the input file regarding BASIS for KIND <"// &
2332 trim(qs_kind%name)//"> will be ignored!")
2333
2334 CASE (do_method_pw)
2335 ! PW DFT
2336 ! Allocate and initialise the potential data set structure
2337 IF (potential_name /= '') THEN
2338 SELECT CASE (trim(potential_type))
2339 CASE ("ALL", "ECP", "CNEO")
2340 CALL cp_abort(__location__, &
2341 "PW DFT calculations only with potential type UPF or GTH possible."// &
2342 " <"//trim(potential_type)//"> was specified "// &
2343 "for the atomic kind <"//trim(qs_kind%name))
2344 CASE ("GTH")
2345 IF (potential_fn_kind == "-") THEN
2346 CALL section_vals_val_get(dft_section, "POTENTIAL_FILE_NAME", c_val=potential_file_name)
2347 ELSE
2348 potential_file_name = potential_fn_kind
2349 END IF
2350 CALL allocate_potential(qs_kind%gth_potential)
2351 CALL read_potential(qs_kind%element_symbol, potential_name, &
2352 qs_kind%gth_potential, zeff_correction, para_env, &
2353 potential_file_name, potential_section, update_input, &
2354 monovalent=qs_kind%monovalent)
2355 CALL set_potential(qs_kind%gth_potential, z=z)
2356 CALL get_qs_kind(qs_kind, elec_conf=elec_conf)
2357 IF (.NOT. ASSOCIATED(elec_conf)) THEN
2358 CALL get_potential(potential=qs_kind%gth_potential, elec_conf=elec_conf)
2359 CALL set_qs_kind(qs_kind, elec_conf=elec_conf)
2360 ELSE
2361 CALL set_potential(potential=qs_kind%gth_potential, elec_conf=elec_conf)
2362 END IF
2363 CASE ("UPF")
2364 ALLOCATE (qs_kind%upf_potential)
2365 qs_kind%upf_potential%zion = 0
2366 qs_kind%upf_potential%filename = adjustl(trim(potential_name))
2367 CALL get_qs_kind(qs_kind, elec_conf=elec_conf)
2368 IF (.NOT. ASSOCIATED(elec_conf)) THEN
2369 CALL set_qs_kind(qs_kind, elec_conf=qs_kind%upf_potential%econf)
2370 END IF
2371 CASE DEFAULT
2372 CALL cp_abort(__location__, &
2373 "An invalid potential type <"// &
2374 trim(potential_type)//"> was specified "// &
2375 "for the atomic kind <"// &
2376 trim(qs_kind%name))
2377 END SELECT
2378 ELSE
2379 CALL cp_abort(__location__, &
2380 "No potential type was defined for the "// &
2381 "atomic kind <"//trim(qs_kind%name)//">")
2382 END IF
2383
2384 CASE DEFAULT
2385
2386 ! set ngauss for STO expansion
2387 CALL section_vals_val_get(dft_section, "QS%STO_NG", i_val=ngauss)
2388 ! Allocate and initialise the basis set data set structure
2389 ! first external basis sets
2390 DO i = 1, nb_rep
2391 SELECT CASE (basis_set_form(i))
2392 CASE ("GTO")
2393 NULLIFY (tmp_basis_set)
2394 CALL allocate_gto_basis_set(tmp_basis_set)
2395 CALL read_gto_basis_set(qs_kind%element_symbol, basis_set_name(i), &
2396 tmp_basis_set, para_env, dft_section)
2397 CASE ("STO")
2398 NULLIFY (sto_basis_set)
2399 CALL allocate_sto_basis_set(sto_basis_set)
2400 CALL read_sto_basis_set(qs_kind%element_symbol, basis_set_name(i), &
2401 sto_basis_set, para_env, dft_section)
2402 NULLIFY (tmp_basis_set)
2403 CALL create_gto_from_sto_basis(sto_basis_set, tmp_basis_set, ngauss)
2404 CALL deallocate_sto_basis_set(sto_basis_set)
2405 CASE DEFAULT
2406 CALL cp_abort(__location__, &
2407 "Invalid basis set form "//trim(basis_set_form(i))// &
2408 "for atomic kind <"//trim(qs_kind%name)//">")
2409 END SELECT
2410 tmp = basis_set_type(i)
2411 CALL uppercase(tmp)
2412 CALL add_basis_set_to_container(qs_kind%basis_sets, tmp_basis_set, tmp)
2413 END DO
2414 ! now explicit basis sets
2415 IF (explicit_basis) THEN
2416 CALL section_vals_get(basis_section, n_repetition=nexp)
2417 DO i = 1, nexp
2418 NULLIFY (tmp_basis_set)
2419 CALL allocate_gto_basis_set(tmp_basis_set)
2420 CALL read_gto_basis_set(qs_kind%element_symbol, basis_type, &
2421 tmp_basis_set, basis_section, i, dft_section)
2422 tmp = basis_type
2423 CALL uppercase(tmp)
2424 CALL add_basis_set_to_container(qs_kind%basis_sets, tmp_basis_set, tmp)
2425 END DO
2426 END IF
2427 ! combine multiple basis sets
2428 DO i = 1, SIZE(qs_kind%basis_sets)
2429 NULLIFY (tmp_basis_set)
2430 CALL get_basis_from_container(qs_kind%basis_sets, basis_set=tmp_basis_set, &
2431 inumbas=i, basis_type=basis_type)
2432 IF (basis_type == "") cycle
2433 jj = i
2434 DO j = i + 1, SIZE(qs_kind%basis_sets)
2435 jj = jj + 1
2436 NULLIFY (sup_basis_set)
2437 CALL get_basis_from_container(qs_kind%basis_sets, basis_set=sup_basis_set, &
2438 inumbas=jj, basis_type=tmp)
2439 IF (basis_type == tmp) THEN
2440 ! we found a match, combine the basis sets and delete the second
2441 CALL combine_basis_sets(tmp_basis_set, sup_basis_set)
2442 CALL remove_basis_from_container(qs_kind%basis_sets, jj)
2443 jj = jj - 1
2444 END IF
2445 END DO
2446 NULLIFY (sup_basis_set)
2447 END DO
2448
2449 ! check that we have an orbital basis set
2450 nobasis = .true.
2451 DO i = 1, SIZE(qs_kind%basis_sets)
2452 NULLIFY (tmp_basis_set)
2453 CALL get_basis_from_container(qs_kind%basis_sets, basis_set=tmp_basis_set, &
2454 inumbas=i, basis_type=basis_type)
2455 IF (basis_type == "ORB") nobasis = .false.
2456 END DO
2457 IF (nobasis) THEN
2458 CALL cp_abort(__location__, &
2459 "No basis set type was defined for the "// &
2460 "atomic kind <"//trim(qs_kind%name)//">")
2461 END IF
2462
2463 ! If Ghost atom we don't need to allocate/initialize anything connected to POTENTIAL
2464 IF (qs_kind%ghost .OR. qs_kind%floating) THEN
2465 IF (ASSOCIATED(qs_kind%elec_conf)) qs_kind%elec_conf = 0
2466 ELSE
2467 ! Allocate and initialise the potential data set structure
2468 IF ((potential_name /= '') .OR. explicit_potential) THEN
2469 ! determine the pseudopotential file to search
2470 IF (potential_fn_kind == "-") THEN
2471 CALL section_vals_val_get(dft_section, "POTENTIAL_FILE_NAME", c_val=potential_file_name)
2472 ELSE
2473 potential_file_name = potential_fn_kind
2474 END IF
2475 !
2476 SELECT CASE (trim(potential_type))
2477 CASE ("ALL")
2478 CALL allocate_potential(qs_kind%all_potential)
2479 CALL read_potential(qs_kind%element_symbol, potential_name, &
2480 qs_kind%all_potential, zeff_correction, para_env, &
2481 potential_file_name, potential_section, update_input)
2482 CALL set_potential(qs_kind%all_potential, z=z)
2483 CALL get_qs_kind(qs_kind, elec_conf=elec_conf)
2484 IF (.NOT. ASSOCIATED(elec_conf)) THEN
2485 CALL get_potential(potential=qs_kind%all_potential, elec_conf=elec_conf)
2486 CALL set_qs_kind(qs_kind, elec_conf=elec_conf)
2487 ELSE
2488 CALL set_potential(potential=qs_kind%all_potential, elec_conf=elec_conf)
2489 END IF
2490 CASE ("GTH")
2491 CALL allocate_potential(qs_kind%gth_potential)
2492 CALL read_potential(qs_kind%element_symbol, potential_name, &
2493 qs_kind%gth_potential, zeff_correction, para_env, &
2494 potential_file_name, potential_section, update_input, &
2495 monovalent=qs_kind%monovalent)
2496 CALL set_potential(qs_kind%gth_potential, z=z)
2497 CALL get_qs_kind(qs_kind, elec_conf=elec_conf)
2498 IF (.NOT. ASSOCIATED(elec_conf)) THEN
2499 CALL get_potential(potential=qs_kind%gth_potential, elec_conf=elec_conf)
2500 CALL set_qs_kind(qs_kind, elec_conf=elec_conf)
2501 ELSE
2502 CALL set_potential(potential=qs_kind%gth_potential, elec_conf=elec_conf)
2503 END IF
2504 CASE ("ECP")
2505 CALL allocate_potential(qs_kind%sgp_potential)
2506 CALL get_potential(qs_kind%sgp_potential, description=description)
2507 CALL read_ecp_potential(ptable(z)%symbol, ecppot, &
2508 potential_name, potential_file_name, potential_section)
2509 IF (ecp_semi_local) THEN
2510 description(1) = "Semi-local Gaussian pseudopotential "
2511 description(2) = "ECP "//trim(potential_name)
2512 description(3) = "LIBGRPP: A. V. Oleynichenko et al., Symmetry 15 197 2023"
2513 description(4) = " "
2514 ELSE
2515 description(4) = "ECP "//trim(potential_name)
2516 END IF
2517 CALL set_potential(qs_kind%sgp_potential, name=ecppot%pname, description=description, &
2518 zeff=ecppot%zion, z=z, ecp_local=.true., ecp_semi_local=ecp_semi_local, &
2519 nloc=ecppot%nloc, nrloc=ecppot%nrloc, aloc=ecppot%aloc, bloc=ecppot%bloc, &
2520 has_nlcc=.false.)
2521 CALL set_potential(qs_kind%sgp_potential, sl_lmax=ecppot%lmax, &
2522 npot=ecppot%npot, nrpot=ecppot%nrpot, apot=ecppot%apot, bpot=ecppot%bpot)
2523 ! convert PP
2524 IF (.NOT. ecp_semi_local) THEN
2525 cpabort("ECPs are only well tested in their semi-local form")
2526 CALL get_qs_kind(qs_kind, basis_set=orb_basis_set)
2527 CALL sgp_construction(sgp_pot=sgppot, ecp_pot=ecppot, orb_basis=orb_basis_set, error=error)
2528 IF (iounit > 0 .AND. .NOT. silent) THEN
2529 WRITE (iounit, "(/,T2,'PP Transformation for ',A)") trim(ecppot%pname)
2530 IF (sgppot%has_local) THEN
2531 WRITE (iounit, "(T8,'Accuracy for local part:',T41,F10.3,'%',T61,F20.12)") error(4), error(1)
2532 END IF
2533 IF (sgppot%has_nonlocal) THEN
2534 WRITE (iounit, "(T8,'Accuracy for nonlocal part:',T41,F10.3,'%',T61,F20.12)") error(5), error(2)
2535 END IF
2536 IF (sgppot%has_nlcc) THEN
2537 WRITE (iounit, "(T8,'Accuracy for NLCC density:',T61,F20.12)") error(3)
2538 END IF
2539 END IF
2540 END IF
2541 IF (sgppot%has_nonlocal) THEN
2542 CALL set_potential(qs_kind%sgp_potential, n_nonlocal=sgppot%n_nonlocal, lmax=sgppot%lmax, &
2543 is_nonlocal=sgppot%is_nonlocal)
2544 nnl = sgppot%n_nonlocal
2545 nppnl = 0
2546 DO l = 0, sgppot%lmax
2547 nppnl = nppnl + nnl*nco(l)
2548 END DO
2549 l = sgppot%lmax
2550 ALLOCATE (a_nl(nnl), h_nl(nnl, 0:l), c_nl(nnl, nnl, 0:l))
2551 a_nl(:) = sgppot%a_nonlocal(:)
2552 h_nl(:, :) = sgppot%h_nonlocal(:, :)
2553 DO l = 0, sgppot%lmax
2554 c_nl(:, :, l) = sgppot%c_nonlocal(:, :, l)*sqrt(2._dp*l + 1.0_dp)
2555 END DO
2556 CALL set_potential(qs_kind%sgp_potential, nppnl=nppnl, a_nonlocal=a_nl, h_nonlocal=h_nl, c_nonlocal=c_nl)
2557 ELSE
2558 CALL set_potential(qs_kind%sgp_potential, n_nonlocal=0, lmax=-1, is_nonlocal=sgppot%is_nonlocal)
2559 CALL set_potential(qs_kind%sgp_potential, nppnl=0)
2560 END IF
2561 !
2562 cpassert(.NOT. sgppot%has_local)
2563 cpassert(.NOT. sgppot%has_nlcc)
2564 ! core
2565 rc = 0.5_dp*qs_kind%covalent_radius*angstrom
2566 rc = max(rc, 0.2_dp)
2567 rc = min(rc, 1.0_dp)
2568 alpha = 1.0_dp/(2.0_dp*rc**2)
2569 ccore = ecppot%zion*sqrt((alpha/pi)**3)
2570 CALL set_potential(qs_kind%sgp_potential, alpha_core_charge=alpha, ccore_charge=ccore, &
2571 core_charge_radius=rc)
2572 CALL atom_sgp_release(sgppot)
2573 CALL get_qs_kind(qs_kind, elec_conf=elec_conf)
2574 IF (.NOT. ASSOCIATED(elec_conf)) THEN
2575 CALL set_qs_kind(qs_kind, elec_conf=ecppot%econf)
2576 END IF
2577 CALL get_qs_kind(qs_kind, elec_conf=elec_conf)
2578 CALL set_potential(qs_kind%sgp_potential, elec_conf=elec_conf)
2579 CASE ("UPF")
2580 CALL allocate_potential(qs_kind%sgp_potential)
2581 CALL get_potential(qs_kind%sgp_potential, description=description)
2582 description(4) = "UPF "//trim(potential_name)
2583 CALL atom_read_upf(upfpot, potential_name)
2584 CALL set_potential(qs_kind%sgp_potential, name=upfpot%pname, description=description, &
2585 zeff=upfpot%zion, z=z, has_nlcc=upfpot%core_correction)
2586 ! convert pp
2587 CALL sgp_construction(sgp_pot=sgppot, upf_pot=upfpot, error=error)
2588 IF (iounit > 0 .AND. .NOT. silent) THEN
2589 WRITE (iounit, "(/,T2,'PP Transformation for ',A)") trim(upfpot%pname)
2590 IF (sgppot%has_local) THEN
2591 WRITE (iounit, "(T8,'Accuracy for local part:',T61,F20.12)") error(1)
2592 END IF
2593 IF (sgppot%has_nonlocal) THEN
2594 WRITE (iounit, "(T8,'Accuracy for nonlocal part:',T61,F20.12)") error(2)
2595 END IF
2596 IF (sgppot%has_nlcc) THEN
2597 WRITE (iounit, "(T8,'Accuracy for NLCC density:',T61,F20.12)") error(3)
2598 END IF
2599 END IF
2600 IF (sgppot%has_nonlocal) THEN
2601 CALL set_potential(qs_kind%sgp_potential, n_nonlocal=sgppot%n_nonlocal, lmax=sgppot%lmax, &
2602 is_nonlocal=sgppot%is_nonlocal)
2603 nnl = sgppot%n_nonlocal
2604 nppnl = 0
2605 DO l = 0, sgppot%lmax
2606 nppnl = nppnl + nnl*nco(l)
2607 END DO
2608 l = sgppot%lmax
2609 ALLOCATE (a_nl(nnl), h_nl(nnl, 0:l), c_nl(nnl, nnl, 0:l))
2610 a_nl(:) = sgppot%a_nonlocal(:)
2611 h_nl(:, :) = sgppot%h_nonlocal(:, :)
2612 c_nl(:, :, :) = sgppot%c_nonlocal(:, :, :)
2613 CALL set_potential(qs_kind%sgp_potential, nppnl=nppnl, a_nonlocal=a_nl, h_nonlocal=h_nl, c_nonlocal=c_nl)
2614 ELSE
2615 CALL set_potential(qs_kind%sgp_potential, n_nonlocal=0, lmax=-1, is_nonlocal=sgppot%is_nonlocal)
2616 CALL set_potential(qs_kind%sgp_potential, nppnl=0)
2617 END IF
2618 cpassert(sgppot%has_local)
2619 ! core
2620 rc = sgppot%ac_local
2621 alpha = 1.0_dp/(2.0_dp*rc**2)
2622 ccore = upfpot%zion*sqrt((alpha/pi)**3)
2623 CALL set_potential(qs_kind%sgp_potential, alpha_core_charge=alpha, ccore_charge=ccore, &
2624 core_charge_radius=rc)
2625 ! local potential
2626 nloc = sgppot%n_local
2627 ALLOCATE (aloc(nloc), cloc(nloc))
2628 aloc(1:nloc) = sgppot%a_local(1:nloc)
2629 cloc(1:nloc) = sgppot%c_local(1:nloc)
2630 CALL set_potential(qs_kind%sgp_potential, n_local=nloc, a_local=aloc, c_local=cloc)
2631 IF (sgppot%has_nlcc) THEN
2632 nlcc = sgppot%n_nlcc
2633 ALLOCATE (anlcc(nlcc), cnlcc(nlcc))
2634 anlcc(1:nlcc) = sgppot%a_nlcc(1:nlcc)
2635 cnlcc(1:nlcc) = sgppot%c_nlcc(1:nlcc)
2636 CALL set_potential(qs_kind%sgp_potential, has_nlcc=.true., n_nlcc=nlcc, a_nlcc=anlcc, c_nlcc=cnlcc)
2637 END IF
2638 CALL set_potential(qs_kind%sgp_potential, z=z)
2639 CALL get_qs_kind(qs_kind, elec_conf=elec_conf)
2640 IF (.NOT. ASSOCIATED(elec_conf)) THEN
2641 CALL set_qs_kind(qs_kind, elec_conf=upfpot%econf)
2642 END IF
2643 CALL get_qs_kind(qs_kind, elec_conf=elec_conf)
2644 CALL set_potential(qs_kind%sgp_potential, elec_conf=elec_conf)
2645 CALL atom_release_upf(upfpot)
2646 CALL atom_sgp_release(sgppot)
2647 CASE ("CNEO")
2648 IF (zeff_correction /= 0.0_dp) &
2649 cpabort("CORE_CORRECTION is not compatible with CNEO")
2650 CALL allocate_cneo_potential(qs_kind%cneo_potential)
2651 CALL set_cneo_potential(qs_kind%cneo_potential, z=z)
2652 mass = 0.0_dp
2653 ! Input mass is the mass of the neutral atom, not the nucleus.
2654 ! The mass of electrons will get subtracted later.
2655 ! In principle, the most abundant pure isotope mass should be used.
2656 IF (k_rep > 0) THEN
2657 CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
2658 keyword_name="MASS", n_rep_val=i)
2659 IF (i > 0) CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
2660 keyword_name="MASS", r_val=mass)
2661 END IF
2662 ! Remove electron mass from atomic mass to get nuclear mass
2663 IF (mass - real(z, dp)*0.000548579909_dp > 0.0_dp) THEN
2664 mass = mass - real(z, dp)*0.000548579909_dp
2665 CALL set_cneo_potential(qs_kind%cneo_potential, mass=mass)
2666 END IF
2667 ! In case the mass is not set by user, get the default mass from z
2668 CALL get_cneo_potential(qs_kind%cneo_potential, mass=mass)
2669 ! Warn if the mass is from ptable
2670 IF (abs(mass + real(z, dp)*0.000548579909_dp - ptable(z)%amass) < 1.e-4_dp) THEN
2671 CALL cp_warn(__location__, &
2672 "Atomic mass of the atomic kind <"//trim(qs_kind%name)// &
2673 "> is very close to its average mass. Is it a pure isotope? "// &
2674 "Pure isotopes are preferable for CNEO. "// &
2675 "(e.g., mass of 1H is 1.007825, not 1.00794)")
2676 END IF
2677 CALL get_qs_kind(qs_kind, elec_conf=elec_conf)
2678 IF (.NOT. ASSOCIATED(elec_conf)) THEN
2679 CALL get_cneo_potential(potential=qs_kind%cneo_potential, elec_conf=elec_conf)
2680 CALL set_qs_kind(qs_kind, elec_conf=elec_conf)
2681 ELSE
2682 CALL set_cneo_potential(potential=qs_kind%cneo_potential, elec_conf=elec_conf)
2683 END IF
2684 CASE DEFAULT
2685 CALL cp_abort(__location__, &
2686 "An invalid potential type <"// &
2687 trim(potential_name)//"> was specified "// &
2688 "for the atomic kind <"// &
2689 trim(qs_kind%name))
2690 END SELECT
2691 ELSE
2692 CALL cp_abort(__location__, &
2693 "No potential type was defined for the "// &
2694 "atomic kind <"//trim(qs_kind%name)//">")
2695 END IF
2696
2697 CALL check_potential_basis_compatibility(qs_kind)
2698
2699 ! Allocate and initialise the potential data set structure
2700 IF ((kgpot_name /= '') .OR. explicit_kgpot) THEN
2701 ipos = index(kgpot_name, "-")
2702 IF (ipos > 1) THEN
2703 kgpot_type = kgpot_name(:ipos - 1)
2704 ELSE
2705 kgpot_type = kgpot_name
2706 END IF
2707 CALL uppercase(kgpot_type)
2708
2709 SELECT CASE (trim(kgpot_type))
2710 CASE ("TNADD")
2711 ! determine the pseudopotential file to search
2712 IF (kg_potential_fn_kind == "-") THEN
2713 CALL section_vals_val_get(dft_section, "POTENTIAL_FILE_NAME", c_val=potential_file_name)
2714 ELSE
2715 potential_file_name = kg_potential_fn_kind
2716 END IF
2717 CALL allocate_potential(qs_kind%tnadd_potential)
2718 CALL read_potential(qs_kind%element_symbol, kgpot_name, &
2719 qs_kind%tnadd_potential, para_env, &
2720 potential_file_name, kgpot_section, update_input)
2721 CASE ("NONE")
2722 NULLIFY (qs_kind%tnadd_potential)
2723 CASE DEFAULT
2724 CALL cp_abort(__location__, &
2725 "An invalid kg_potential type <"// &
2726 trim(potential_name)//"> was specified "// &
2727 "for the atomic kind <"// &
2728 trim(qs_kind%name))
2729 END SELECT
2730 END IF
2731 END IF
2732
2733 ! check that we have a nuclear orbital basis set if CNEO is requested
2734 nobasis_nuc = ASSOCIATED(qs_kind%cneo_potential)
2735 DO i = 1, SIZE(qs_kind%basis_sets)
2736 NULLIFY (tmp_basis_set)
2737 CALL get_basis_from_container(qs_kind%basis_sets, basis_set=tmp_basis_set, &
2738 inumbas=i, basis_type=basis_type)
2739 IF (basis_type == "NUC") THEN
2740 nobasis_nuc = .false.
2741 IF (.NOT. ASSOCIATED(qs_kind%cneo_potential)) THEN
2742 CALL cp_warn(__location__, &
2743 "POTENTIAL is not set to CNEO, NUC type basis set for KIND <"// &
2744 trim(qs_kind%name)//"> will be ignored!")
2745 END IF
2746 END IF
2747 END DO
2748 IF (nobasis_nuc) THEN
2749 CALL cp_abort(__location__, &
2750 "No NUC type basis set was defined for the "// &
2751 "atomic kind <"//trim(qs_kind%name)//">, which is required by "// &
2752 "POTENTIAL CNEO.")
2753 END IF
2754 END SELECT
2755
2756 CALL timestop(handle)
2757
2758 END SUBROUTINE read_qs_kind
2759
2760! **************************************************************************************************
2761!> \brief Ensure pseudo-potential and basis set were optimized for same number of valence electrons
2762!> \param qs_kind ...
2763!> \author Ole Schuett
2764! **************************************************************************************************
2765 SUBROUTINE check_potential_basis_compatibility(qs_kind)
2766 TYPE(qs_kind_type), INTENT(INOUT) :: qs_kind
2767
2768 CHARACTER(LEN=default_string_length) :: name
2769 INTEGER :: nbs, npp
2770 TYPE(gth_potential_type), POINTER :: gth_potential
2771 TYPE(gto_basis_set_type), POINTER :: basis_set
2772
2773 CALL get_qs_kind(qs_kind, name=name, gth_potential=gth_potential, basis_set=basis_set)
2774
2775 npp = -1; nbs = -1
2776 IF (ASSOCIATED(gth_potential)) &
2777 npp = parse_valence_electrons(gth_potential%aliases)
2778 IF (ASSOCIATED(basis_set)) &
2779 nbs = parse_valence_electrons(basis_set%aliases)
2780
2781 IF (npp >= 0 .AND. nbs >= 0 .AND. npp /= nbs) &
2782 CALL cp_abort(__location__, "Basis-set and pseudo-potential of atomic kind '"//trim(name)//"'"// &
2783 " were optimized for different valence electron numbers.")
2784
2785 END SUBROUTINE check_potential_basis_compatibility
2786
2787! **************************************************************************************************
2788!> \brief Tries to parse valence eletron number using "-QXXX" notation, returns -1 if not found.
2789!> \param string ...
2790!> \return ...
2791!> \author Ole Schuett
2792! **************************************************************************************************
2793 FUNCTION parse_valence_electrons(string) RESULT(n)
2794 CHARACTER(*) :: string
2795 INTEGER :: n
2796
2797 INTEGER :: i, istat, j
2798
2799 i = index(string, "-Q", .true.)
2800 IF (i == 0) THEN
2801 n = -1
2802 ELSE
2803 j = scan(string(i + 2:), "- ")
2804 READ (string(i + 2:i + j), '(I3)', iostat=istat) n
2805 IF (istat /= 0) n = -1
2806 END IF
2807
2808 END FUNCTION parse_valence_electrons
2809
2810! **************************************************************************************************
2811!> \brief Read an atomic kind set data set from the input file.
2812!> \param qs_kind_set ...
2813!> \param atomic_kind_set ...
2814!> \param kind_section ...
2815!> \param para_env ...
2816!> \param force_env_section ...
2817!> \param silent ...
2818! **************************************************************************************************
2819 SUBROUTINE create_qs_kind_set(qs_kind_set, atomic_kind_set, kind_section, para_env, &
2820 force_env_section, silent)
2821
2822 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
2823 TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
2824 TYPE(section_vals_type), POINTER :: kind_section
2825 TYPE(mp_para_env_type), POINTER :: para_env
2826 TYPE(section_vals_type), POINTER :: force_env_section
2827 LOGICAL, INTENT(IN) :: silent
2828
2829 CHARACTER(len=*), PARAMETER :: routinen = 'create_qs_kind_set'
2830
2831 INTEGER :: handle, ikind, method, nkind, qs_method
2832 LOGICAL :: no_fail
2833
2834 CALL timeset(routinen, handle)
2835
2836 IF (ASSOCIATED(qs_kind_set)) cpabort("create_qs_kind_set: qs_kind_set already associated")
2837 IF (.NOT. ASSOCIATED(atomic_kind_set)) cpabort("create_qs_kind_set: atomic_kind_set not associated")
2838
2839 no_fail = .false.
2840
2841 ! Between all methods only SE and DFTB/xTB may not need a KIND section.
2842 CALL section_vals_val_get(force_env_section, "METHOD", i_val=method)
2843 IF (method == do_qs) THEN
2844 CALL section_vals_val_get(force_env_section, "DFT%QS%METHOD", i_val=qs_method)
2845 SELECT CASE (qs_method)
2848 no_fail = .true.
2849 CASE (do_method_dftb)
2850 no_fail = .true.
2851 CASE (do_method_xtb)
2852 no_fail = .true.
2853 END SELECT
2854 ELSE IF (method == do_sirius) THEN
2855 qs_method = do_method_pw
2856 ELSE
2857 qs_method = method
2858 END IF
2859
2860 nkind = SIZE(atomic_kind_set)
2861 ALLOCATE (qs_kind_set(nkind))
2862
2863 DO ikind = 1, nkind
2864 qs_kind_set(ikind)%name = atomic_kind_set(ikind)%name
2865 qs_kind_set(ikind)%element_symbol = atomic_kind_set(ikind)%element_symbol
2866 qs_kind_set(ikind)%natom = atomic_kind_set(ikind)%natom
2867 CALL read_qs_kind(qs_kind_set(ikind), kind_section, para_env, force_env_section, &
2868 no_fail, qs_method, silent)
2869 END DO
2870
2871 CALL timestop(handle)
2872
2873 END SUBROUTINE create_qs_kind_set
2874
2875! **************************************************************************************************
2876!> \brief This routines should perform only checks. no settings are allowed at
2877!> this level anymore..
2878!> \param qs_kind ...
2879!> \param dft_control ...
2880!> \param subsys_section ...
2881! **************************************************************************************************
2882 SUBROUTINE check_qs_kind(qs_kind, dft_control, subsys_section)
2883
2884 TYPE(qs_kind_type), POINTER :: qs_kind
2885 TYPE(dft_control_type), INTENT(IN) :: dft_control
2886 TYPE(section_vals_type), POINTER :: subsys_section
2887
2888 INTEGER :: gfn_type
2889 LOGICAL :: defined
2890 TYPE(qs_dftb_atom_type), POINTER :: dftb_parameter
2891 TYPE(semi_empirical_type), POINTER :: se_parameter
2892 TYPE(xtb_atom_type), POINTER :: xtb_parameter
2893
2894 IF (dft_control%qs_control%semi_empirical) THEN
2895 CALL get_qs_kind(qs_kind, se_parameter=se_parameter)
2896 cpassert(ASSOCIATED(se_parameter))
2897 CALL get_se_param(se_parameter, defined=defined)
2898 cpassert(defined)
2899 CALL write_se_param(se_parameter, subsys_section)
2900 ELSE IF (dft_control%qs_control%dftb) THEN
2901 CALL get_qs_kind(qs_kind, dftb_parameter=dftb_parameter)
2902 cpassert(ASSOCIATED(dftb_parameter))
2903 CALL get_dftb_atom_param(dftb_parameter, defined=defined)
2904 cpassert(defined)
2905 CALL write_dftb_atom_param(dftb_parameter, subsys_section)
2906 ELSE IF (dft_control%qs_control%xtb) THEN
2907 IF (.NOT. (dft_control%qs_control%xtb_control%do_tblite)) THEN
2908 CALL get_qs_kind(qs_kind, xtb_parameter=xtb_parameter)
2909 cpassert(ASSOCIATED(xtb_parameter))
2910 gfn_type = dft_control%qs_control%xtb_control%gfn_type
2911 CALL write_xtb_atom_param(xtb_parameter, gfn_type, subsys_section)
2912 END IF
2913 END IF
2914
2915 END SUBROUTINE check_qs_kind
2916
2917! **************************************************************************************************
2918!> \brief ...
2919!> \param qs_kind_set ...
2920!> \param dft_control ...
2921!> \param subsys_section ...
2922! **************************************************************************************************
2923 SUBROUTINE check_qs_kind_set(qs_kind_set, dft_control, subsys_section)
2924
2925 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
2926 TYPE(dft_control_type), INTENT(IN) :: dft_control
2927 TYPE(section_vals_type), POINTER :: subsys_section
2928
2929 CHARACTER(len=*), PARAMETER :: routinen = 'check_qs_kind_set'
2930
2931 INTEGER :: handle, ikind, nkind
2932 TYPE(qs_kind_type), POINTER :: qs_kind
2933
2934 CALL timeset(routinen, handle)
2935 IF (ASSOCIATED(qs_kind_set)) THEN
2936 nkind = SIZE(qs_kind_set)
2937 DO ikind = 1, nkind
2938 qs_kind => qs_kind_set(ikind)
2939 CALL check_qs_kind(qs_kind, dft_control, subsys_section)
2940 END DO
2941 IF (dft_control%qs_control%xtb) THEN
2942 CALL write_xtb_kab_param(qs_kind_set, subsys_section, &
2943 dft_control%qs_control%xtb_control)
2944 END IF
2945 ELSE
2946 cpabort("The pointer qs_kind_set is not associated")
2947 END IF
2948 CALL timestop(handle)
2949 END SUBROUTINE check_qs_kind_set
2950
2951! **************************************************************************************************
2952!> \brief ...
2953!> \param qs_kind_set ...
2954!> \param subsys_section ...
2955!> \param xtb_control ...
2956! **************************************************************************************************
2957 SUBROUTINE write_xtb_kab_param(qs_kind_set, subsys_section, xtb_control)
2958
2959 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
2960 TYPE(section_vals_type), POINTER :: subsys_section
2961 TYPE(xtb_control_type), POINTER :: xtb_control
2962
2963 CHARACTER(LEN=default_string_length) :: aname, bname
2964 INTEGER :: ikind, io_unit, jkind, nkind, za, zb
2965 TYPE(cp_logger_type), POINTER :: logger
2966 TYPE(qs_kind_type), POINTER :: qs_kinda, qs_kindb
2967 TYPE(xtb_atom_type), POINTER :: xtb_parameter_a, xtb_parameter_b
2968
2969 NULLIFY (logger)
2970 logger => cp_get_default_logger()
2971 IF (btest(cp_print_key_should_output(logger%iter_info, subsys_section, &
2972 "PRINT%KINDS/POTENTIAL"), cp_p_file)) THEN
2973
2974 io_unit = cp_print_key_unit_nr(logger, subsys_section, "PRINT%KINDS", extension=".Log")
2975 IF (io_unit > 0) THEN
2976
2977 WRITE (io_unit, "(/,T2,A)") "xTB| Kab parameters"
2978 nkind = SIZE(qs_kind_set)
2979 DO ikind = 1, nkind
2980 qs_kinda => qs_kind_set(ikind)
2981 CALL get_qs_kind(qs_kinda, xtb_parameter=xtb_parameter_a)
2982 CALL get_xtb_atom_param(xtb_parameter_a, aname=aname, z=za)
2983 DO jkind = ikind, nkind
2984 qs_kindb => qs_kind_set(jkind)
2985 CALL get_qs_kind(qs_kindb, xtb_parameter=xtb_parameter_b)
2986 CALL get_xtb_atom_param(xtb_parameter_b, aname=bname, z=zb)
2987 WRITE (io_unit, "(A,T10,A15,T25,A15,T71,F10.3)") &
2988 " Kab:", trim(aname), trim(bname), xtb_set_kab(za, zb, xtb_control)
2989 END DO
2990 END DO
2991 WRITE (io_unit, *)
2992
2993 END IF
2994
2995 CALL cp_print_key_finished_output(io_unit, logger, subsys_section, "PRINT%KINDS")
2996 END IF
2997
2998 END SUBROUTINE write_xtb_kab_param
2999
3000! **************************************************************************************************
3001!> \brief Set the components of an atomic kind data set.
3002!> \param qs_kind ...
3003!> \param paw_atom ...
3004!> \param ghost ...
3005!> \param floating ...
3006!> \param hard_radius ...
3007!> \param hard0_radius ...
3008!> \param covalent_radius ...
3009!> \param vdw_radius ...
3010!> \param lmax_rho0 ...
3011!> \param zeff ...
3012!> \param no_optimize ...
3013!> \param dispersion ...
3014!> \param u_minus_j ...
3015!> \param reltmat ...
3016!> \param dftb_parameter ...
3017!> \param xtb_parameter ...
3018!> \param elec_conf ...
3019!> \param pao_basis_size ...
3020! **************************************************************************************************
3021 SUBROUTINE set_qs_kind(qs_kind, paw_atom, ghost, floating, hard_radius, hard0_radius, &
3022 covalent_radius, vdw_radius, lmax_rho0, zeff, &
3023 no_optimize, dispersion, u_minus_j, reltmat, &
3024 dftb_parameter, xtb_parameter, &
3025 elec_conf, pao_basis_size)
3026
3027 TYPE(qs_kind_type), INTENT(INOUT) :: qs_kind
3028 LOGICAL, INTENT(IN), OPTIONAL :: paw_atom, ghost, floating
3029 REAL(kind=dp), INTENT(IN), OPTIONAL :: hard_radius, hard0_radius, &
3030 covalent_radius, vdw_radius
3031 INTEGER, INTENT(IN), OPTIONAL :: lmax_rho0
3032 REAL(kind=dp), INTENT(IN), OPTIONAL :: zeff
3033 LOGICAL, INTENT(IN), OPTIONAL :: no_optimize
3034 TYPE(qs_atom_dispersion_type), OPTIONAL, POINTER :: dispersion
3035 REAL(kind=dp), INTENT(IN), OPTIONAL :: u_minus_j
3036 REAL(kind=dp), DIMENSION(:, :), OPTIONAL, POINTER :: reltmat
3037 TYPE(qs_dftb_atom_type), OPTIONAL, POINTER :: dftb_parameter
3038 TYPE(xtb_atom_type), OPTIONAL, POINTER :: xtb_parameter
3039 INTEGER, DIMENSION(:), INTENT(IN), OPTIONAL :: elec_conf
3040 INTEGER, INTENT(IN), OPTIONAL :: pao_basis_size
3041
3042 IF (PRESENT(dftb_parameter)) qs_kind%dftb_parameter => dftb_parameter
3043 IF (PRESENT(xtb_parameter)) qs_kind%xtb_parameter => xtb_parameter
3044 IF (PRESENT(elec_conf)) THEN
3045 IF (ASSOCIATED(qs_kind%elec_conf)) THEN
3046 DEALLOCATE (qs_kind%elec_conf)
3047 END IF
3048 ALLOCATE (qs_kind%elec_conf(0:SIZE(elec_conf) - 1))
3049 qs_kind%elec_conf(:) = elec_conf(:)
3050 END IF
3051 IF (PRESENT(paw_atom)) qs_kind%paw_atom = paw_atom
3052 IF (PRESENT(hard_radius)) qs_kind%hard_radius = hard_radius
3053 IF (PRESENT(hard0_radius)) qs_kind%hard0_radius = hard0_radius
3054 IF (PRESENT(covalent_radius)) qs_kind%covalent_radius = covalent_radius
3055 IF (PRESENT(vdw_radius)) qs_kind%vdw_radius = vdw_radius
3056 IF (PRESENT(lmax_rho0)) qs_kind%lmax_rho0 = lmax_rho0
3057 IF (PRESENT(zeff)) THEN
3058 IF (ASSOCIATED(qs_kind%all_potential)) THEN
3059 CALL set_potential(potential=qs_kind%all_potential, zeff=zeff)
3060 ELSE IF (ASSOCIATED(qs_kind%gth_potential)) THEN
3061 CALL set_potential(potential=qs_kind%gth_potential, zeff=zeff)
3062 ELSE IF (ASSOCIATED(qs_kind%sgp_potential)) THEN
3063 CALL set_potential(potential=qs_kind%sgp_potential, zeff=zeff)
3064 ELSE IF (ASSOCIATED(qs_kind%cneo_potential)) THEN
3065 cpabort("CNEO potential ZEFF should not be manually set")
3066 END IF
3067 END IF
3068 IF (PRESENT(ghost)) qs_kind%ghost = ghost
3069
3070 IF (PRESENT(floating)) qs_kind%floating = floating
3071
3072 IF (PRESENT(no_optimize)) qs_kind%no_optimize = no_optimize
3073
3074 IF (PRESENT(dispersion)) qs_kind%dispersion => dispersion
3075
3076 IF (PRESENT(u_minus_j)) THEN
3077 IF (ASSOCIATED(qs_kind%dft_plus_u)) THEN
3078 qs_kind%dft_plus_u%u_minus_j = u_minus_j
3079 END IF
3080 END IF
3081
3082 IF (PRESENT(reltmat)) qs_kind%reltmat => reltmat
3083
3084 IF (PRESENT(pao_basis_size)) qs_kind%pao_basis_size = pao_basis_size
3085
3086 END SUBROUTINE set_qs_kind
3087
3088! **************************************************************************************************
3089!> \brief Write an atomic kind data set to the output unit.
3090!> \param qs_kind ...
3091!> \param kind_number ...
3092!> \param output_unit ...
3093!> \par History
3094!> Creation (09.02.2002,MK)
3095! **************************************************************************************************
3096 SUBROUTINE write_qs_kind(qs_kind, kind_number, output_unit)
3097
3098 TYPE(qs_kind_type), POINTER :: qs_kind
3099 INTEGER, INTENT(in) :: kind_number, output_unit
3100
3101 CHARACTER(LEN=3) :: yon
3102 CHARACTER(LEN=default_string_length) :: basis_type, bstring
3103 INTEGER :: ibas
3104 LOGICAL :: do_print
3105 TYPE(gto_basis_set_type), POINTER :: tmp_basis
3106
3107 IF (output_unit > 0) THEN
3108
3109 IF (ASSOCIATED(qs_kind)) THEN
3110 WRITE (unit=output_unit, fmt="(/,T2,I2,A,T57,A,T75,I6)") &
3111 kind_number, ". Atomic kind: "//trim(qs_kind%name), &
3112 "Number of atoms: ", qs_kind%natom
3113
3114 DO ibas = 1, SIZE(qs_kind%basis_sets, 1)
3115 NULLIFY (tmp_basis)
3116 CALL get_basis_from_container(qs_kind%basis_sets, basis_set=tmp_basis, &
3117 inumbas=ibas, basis_type=basis_type)
3118 do_print = .true.
3119 SELECT CASE (basis_type)
3120 CASE DEFAULT
3121 bstring = "Basis Set"
3122 do_print = .false.
3123 CASE ("ORB")
3124 bstring = "Orbital Basis Set"
3125 CASE ("ORB_SOFT")
3126 bstring = "GAPW Soft Basis Set"
3127 do_print = .false.
3128 CASE ("AUX")
3129 bstring = "Auxiliary Basis Set"
3130 CASE ("MIN")
3131 bstring = "Minimal Basis Set"
3132 CASE ("RI_AUX")
3133 bstring = "RI Auxiliary Basis Set"
3134 CASE ("AUX_FIT")
3135 bstring = "Auxiliary Fit Basis Set"
3136 CASE ("LRI_AUX")
3137 bstring = "LRI Basis Set"
3138 CASE ("P_LRI_AUX")
3139 bstring = "LRI Basis Set for TDDFPT"
3140 CASE ("RI_XAS")
3141 bstring = "RI XAS Basis Set"
3142 CASE ("RI_HFX")
3143 bstring = "RI HFX Basis Set"
3144 CASE ("NUC")
3145 bstring = "Nuclear Basis Set"
3146 do_print = .false.
3147 CASE ("NUC_SOFT")
3148 bstring = "Nuclear Soft Basis Set"
3149 do_print = .false.
3150 END SELECT
3151
3152 IF (do_print) THEN
3153 CALL write_orb_basis_set(tmp_basis, output_unit, bstring)
3154 END IF
3155
3156 END DO
3157
3158 IF (qs_kind%ghost) THEN
3159 WRITE (unit=output_unit, fmt="(/,T6,A)") &
3160 "The atoms of this atomic kind are GHOST atoms!"
3161 END IF
3162 IF (qs_kind%monovalent) THEN
3163 WRITE (unit=output_unit, fmt="(/,T6,A)") &
3164 "The atoms of this atomic kind are MONOVALENT!"
3165 END IF
3166 IF (qs_kind%floating) THEN
3167 WRITE (unit=output_unit, fmt="(/,T6,A)") &
3168 "The atoms of this atomic kind are FLOATING BASIS FUNCTIONS."
3169 END IF
3170 IF (qs_kind%covalent_radius > 0.0_dp) THEN
3171 WRITE (unit=output_unit, fmt="(/,T8,A,T71,F10.3)") &
3172 "Atomic covalent radius [Angstrom]:", &
3173 qs_kind%covalent_radius*angstrom
3174 END IF
3175 IF (qs_kind%vdw_radius > 0.0_dp) THEN
3176 WRITE (unit=output_unit, fmt="(/,T8,A,T71,F10.3)") &
3177 "Atomic van der Waals radius [Angstrom]:", &
3178 qs_kind%vdw_radius*angstrom
3179 END IF
3180 IF (qs_kind%paw_atom) THEN
3181 WRITE (unit=output_unit, fmt="(/,T6,A)") &
3182 "The atoms of this atomic kind are PAW atoms (GAPW):"
3183 WRITE (unit=output_unit, fmt="(T8,A,T71,F10.3)") &
3184 "Hard Gaussian function radius:", qs_kind%hard_radius, &
3185 "Rho0 radius:", qs_kind%hard0_radius, &
3186 "Maximum GTO radius used for PAW projector construction:", &
3187 qs_kind%max_rad_local
3188 NULLIFY (tmp_basis)
3189 CALL get_basis_from_container(qs_kind%basis_sets, basis_set=tmp_basis, &
3190 basis_type="ORB_SOFT")
3191 CALL write_orb_basis_set(tmp_basis, output_unit, "GAPW Soft Basis Set")
3192 END IF
3193 ! Potentials
3194 IF (ASSOCIATED(qs_kind%all_potential)) CALL write_potential(qs_kind%all_potential, output_unit)
3195 IF (ASSOCIATED(qs_kind%gth_potential)) CALL write_potential(qs_kind%gth_potential, output_unit)
3196 IF (ASSOCIATED(qs_kind%sgp_potential)) CALL write_potential(qs_kind%sgp_potential, output_unit)
3197 IF (ASSOCIATED(qs_kind%tnadd_potential)) CALL write_potential(qs_kind%tnadd_potential, output_unit)
3198 IF (ASSOCIATED(qs_kind%dft_plus_u)) THEN
3199 WRITE (unit=output_unit, fmt="(/,T6,A,/,T8,A,T76,I5,/,T8,A,T73,F8.3)") &
3200 "A DFT+U correction is applied to atoms of this atomic kind:", &
3201 "Angular quantum momentum number L:", qs_kind%dft_plus_u%l, &
3202 "U(eff) = (U - J) value in [eV]:", qs_kind%dft_plus_u%u_minus_j_target*evolt
3203 IF (qs_kind%dft_plus_u%u_ramping > 0.0_dp) THEN
3204 IF (qs_kind%dft_plus_u%init_u_ramping_each_scf) THEN
3205 yon = "YES"
3206 ELSE
3207 yon = " NO"
3208 END IF
3209 WRITE (unit=output_unit, fmt="(T8,A,T73,F8.3,/,T8,A,T73,ES8.1,/,T8,A,T78,A3)") &
3210 "Increment for U ramping in [eV]:", qs_kind%dft_plus_u%u_ramping*evolt, &
3211 "SCF threshold value for U ramping:", qs_kind%dft_plus_u%eps_u_ramping, &
3212 "Set U ramping value to zero before each wavefunction optimisation:", yon
3213 END IF
3214 IF (ASSOCIATED(qs_kind%dft_plus_u%orbitals)) THEN
3215 WRITE (unit=output_unit, fmt="(T8,A)") &
3216 "An initial orbital occupation is requested:"
3217 IF (ASSOCIATED(qs_kind%dft_plus_u%nelec)) THEN
3218 IF (any(qs_kind%dft_plus_u%nelec(:) >= 0.5_dp)) THEN
3219 IF (SIZE(qs_kind%dft_plus_u%nelec) > 1) THEN
3220 WRITE (unit=output_unit, fmt="(T9,A,T75,F6.2)") &
3221 "Number of alpha electrons:", &
3222 qs_kind%dft_plus_u%nelec(1), &
3223 "Number of beta electrons:", &
3224 qs_kind%dft_plus_u%nelec(2)
3225 ELSE
3226 WRITE (unit=output_unit, fmt="(T9,A,T75,F6.2)") &
3227 "Number of electrons:", &
3228 qs_kind%dft_plus_u%nelec(1)
3229 END IF
3230 END IF
3231 END IF
3232 WRITE (unit=output_unit, fmt="(T9,A,(T78,I3))") &
3233 "Preferred (initial) orbital occupation order (orbital M values):", &
3234 qs_kind%dft_plus_u%orbitals(:)
3235 WRITE (unit=output_unit, fmt="(T9,A,T71,ES10.3,/,T9,A,T76,I5)") &
3236 "Threshold value for the SCF convergence criterion:", &
3237 qs_kind%dft_plus_u%eps_scf, &
3238 "Number of initial SCF iterations:", &
3239 qs_kind%dft_plus_u%max_scf
3240 IF (qs_kind%dft_plus_u%smear) THEN
3241 WRITE (unit=output_unit, fmt="(T9,A)") &
3242 "A smearing of the orbital occupations will be performed"
3243 END IF
3244 END IF
3245 END IF
3246 IF (ASSOCIATED(qs_kind%cneo_potential)) THEN
3247 WRITE (unit=output_unit, fmt="(/,T6,A)") &
3248 "The nuclei of this atomic kind are quantum mechanical (CNEO)"
3249 CALL write_cneo_potential(qs_kind%cneo_potential, output_unit)
3250 NULLIFY (tmp_basis)
3251 CALL get_basis_from_container(qs_kind%basis_sets, basis_set=tmp_basis, &
3252 basis_type="NUC")
3253 CALL write_orb_basis_set(tmp_basis, output_unit, "Nuclear Basis Set")
3254 NULLIFY (tmp_basis)
3255 CALL get_basis_from_container(qs_kind%basis_sets, basis_set=tmp_basis, &
3256 basis_type="NUC_SOFT")
3257 CALL write_orb_basis_set(tmp_basis, output_unit, "Nuclear Soft Basis Set")
3258 END IF
3259 ELSE
3260 cpabort("")
3261 END IF
3262
3263 END IF
3264
3265 END SUBROUTINE write_qs_kind
3266
3267! **************************************************************************************************
3268!> \brief Write an atomic kind set data set to the output unit.
3269!> \param qs_kind_set ...
3270!> \param subsys_section ...
3271!> \par History
3272!> Creation (09.02.2002,MK)
3273! **************************************************************************************************
3274 SUBROUTINE write_qs_kind_set(qs_kind_set, subsys_section)
3275 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
3276 TYPE(section_vals_type), POINTER :: subsys_section
3277
3278 CHARACTER(len=*), PARAMETER :: routinen = 'write_qs_kind_set'
3279
3280 INTEGER :: handle, ikind, nkind, output_unit
3281 TYPE(cp_logger_type), POINTER :: logger
3282 TYPE(qs_kind_type), POINTER :: qs_kind
3283
3284 CALL timeset(routinen, handle)
3285
3286 NULLIFY (logger)
3287 logger => cp_get_default_logger()
3288 output_unit = cp_print_key_unit_nr(logger, subsys_section, &
3289 "PRINT%KINDS", extension=".Log")
3290 IF (output_unit > 0) THEN
3291 IF (ASSOCIATED(qs_kind_set)) THEN
3292 WRITE (unit=output_unit, fmt="(/,/,T2,A)") "ATOMIC KIND INFORMATION"
3293 nkind = SIZE(qs_kind_set)
3294 DO ikind = 1, nkind
3295 qs_kind => qs_kind_set(ikind)
3296 CALL write_qs_kind(qs_kind, ikind, output_unit)
3297 END DO
3298 ELSE
3299 cpabort("")
3300 END IF
3301 END IF
3302
3303 CALL cp_print_key_finished_output(output_unit, logger, subsys_section, &
3304 "PRINT%KINDS")
3305
3306 CALL timestop(handle)
3307
3308 END SUBROUTINE write_qs_kind_set
3309
3310! **************************************************************************************************
3311!> \brief Write all the GTO basis sets of an atomic kind set to the output
3312!> unit (for the printing of the unnormalized basis sets as read from
3313!> database).
3314!> \param qs_kind_set ...
3315!> \param subsys_section ...
3316!> \par History
3317!> Creation (17.01.2002,MK)
3318! **************************************************************************************************
3319 SUBROUTINE write_gto_basis_sets(qs_kind_set, subsys_section)
3320
3321 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
3322 TYPE(section_vals_type), POINTER :: subsys_section
3323
3324 CHARACTER(LEN=*), PARAMETER :: routinen = 'write_gto_basis_sets'
3325
3326 CHARACTER(LEN=default_string_length) :: basis_type, bstring
3327 INTEGER :: handle, ibas, ikind, nkind, output_unit
3328 TYPE(cp_logger_type), POINTER :: logger
3329 TYPE(gto_basis_set_type), POINTER :: tmp_basis
3330 TYPE(qs_kind_type), POINTER :: qs_kind
3331
3332 CALL timeset(routinen, handle)
3333
3334 NULLIFY (logger)
3335 logger => cp_get_default_logger()
3336 output_unit = cp_print_key_unit_nr(logger, subsys_section, &
3337 "PRINT%KINDS/BASIS_SET", &
3338 extension=".Log")
3339 IF (output_unit > 0) THEN
3340 IF (ASSOCIATED(qs_kind_set)) THEN
3341 WRITE (unit=output_unit, fmt="(/,/,T2,A)") &
3342 "BASIS SET INFORMATION (Unnormalised Gaussian-type functions)"
3343 nkind = SIZE(qs_kind_set)
3344 DO ikind = 1, nkind
3345 qs_kind => qs_kind_set(ikind)
3346 WRITE (unit=output_unit, fmt="(/,T2,I2,A)") &
3347 ikind, ". Atomic kind: "//trim(qs_kind%name)
3348
3349 DO ibas = 1, SIZE(qs_kind%basis_sets, 1)
3350 NULLIFY (tmp_basis)
3351 CALL get_basis_from_container(qs_kind%basis_sets, basis_set=tmp_basis, &
3352 inumbas=ibas, basis_type=basis_type)
3353 IF (basis_type == "") cycle
3354 SELECT CASE (basis_type)
3355 CASE DEFAULT
3356 bstring = "Basis Set"
3357 CASE ("ORB")
3358 bstring = "Orbital Basis Set"
3359 CASE ("ORB_SOFT")
3360 bstring = "GAPW Soft Basis Set"
3361 CASE ("AUX")
3362 bstring = "Auxiliary Basis Set"
3363 CASE ("MIN")
3364 bstring = "Minimal Basis Set"
3365 CASE ("RI_AUX")
3366 bstring = "RI Auxiliary Basis Set"
3367 CASE ("AUX_FIT")
3368 bstring = "Auxiliary Fit Basis Set"
3369 CASE ("LRI_AUX")
3370 bstring = "LRI Basis Set"
3371 CASE ("P_LRI_AUX")
3372 bstring = "LRI Basis Set for TDDFPT"
3373 CASE ("RI_HFX")
3374 bstring = "RI HFX Basis Set"
3375 CASE ("NUC")
3376 bstring = "Nuclear Basis Set"
3377 IF (.NOT. ASSOCIATED(qs_kind%cneo_potential)) NULLIFY (tmp_basis)
3378 CASE ("NUC_SOFT")
3379 bstring = "Nuclear Soft Basis Set"
3380 IF (.NOT. ASSOCIATED(qs_kind%cneo_potential)) NULLIFY (tmp_basis)
3381 END SELECT
3382
3383 IF (ASSOCIATED(tmp_basis)) CALL write_gto_basis_set(tmp_basis, output_unit, bstring)
3384
3385 END DO
3386
3387 END DO
3388 ELSE
3389 cpabort("")
3390 END IF
3391 END IF
3392
3393 CALL cp_print_key_finished_output(output_unit, logger, subsys_section, &
3394 "PRINT%KINDS/BASIS_SET")
3395
3396 CALL timestop(handle)
3397
3398 END SUBROUTINE write_gto_basis_sets
3399
3400! **************************************************************************************************
3401!> \brief ...
3402!> \param atomic_kind ...
3403!> \param qs_kind ...
3404!> \param ncalc ...
3405!> \param ncore ...
3406!> \param nelem ...
3407!> \param edelta ...
3408! **************************************************************************************************
3409 SUBROUTINE init_atom_electronic_state(atomic_kind, qs_kind, ncalc, ncore, nelem, edelta)
3410
3411 TYPE(atomic_kind_type), INTENT(IN) :: atomic_kind
3412 TYPE(qs_kind_type), INTENT(IN) :: qs_kind
3413 INTEGER, DIMENSION(0:lmat, 10), INTENT(OUT) :: ncalc, ncore, nelem
3414 REAL(kind=dp), DIMENSION(0:lmat, 10, 2), &
3415 INTENT(OUT) :: edelta
3416
3417 INTEGER :: i, ii, is, l, ll, ne, nn, z
3418 INTEGER, DIMENSION(:), POINTER :: econf
3419 INTEGER, DIMENSION(:, :), POINTER :: addel, laddel, naddel
3420 LOGICAL :: bs_occupation, monovalent
3421 REAL(kind=dp) :: dmag, magnetization
3422 TYPE(gth_potential_type), POINTER :: gth_potential
3423 TYPE(sgp_potential_type), POINTER :: sgp_potential
3424
3425 CALL get_atomic_kind(atomic_kind, z=z)
3426 NULLIFY (gth_potential)
3427 CALL get_qs_kind(qs_kind, &
3428 gth_potential=gth_potential, &
3429 sgp_potential=sgp_potential, &
3430 magnetization=magnetization, &
3431 bs_occupation=bs_occupation, &
3432 monovalent=monovalent, &
3433 addel=addel, laddel=laddel, naddel=naddel)
3434
3435 ! electronic state
3436 nelem = 0
3437 ncore = 0
3438 ncalc = 0
3439 edelta = 0.0_dp
3440 IF (monovalent) THEN
3441 ncalc(0, 1) = 1
3442 nelem(0, 1) = 1
3443 ELSE IF (ASSOCIATED(gth_potential)) THEN
3444 CALL get_potential(gth_potential, elec_conf=econf)
3445 CALL set_pseudo_state(econf, z, ncalc, ncore, nelem)
3446 ELSE IF (ASSOCIATED(sgp_potential)) THEN
3447 CALL get_potential(sgp_potential, elec_conf=econf)
3448 CALL set_pseudo_state(econf, z, ncalc, ncore, nelem)
3449 ELSE
3450 DO l = 0, min(lmat, ubound(ptable(z)%e_conv, 1))
3451 ll = 2*(2*l + 1)
3452 nn = ptable(z)%e_conv(l)
3453 ii = 0
3454 DO
3455 ii = ii + 1
3456 IF (nn <= ll) THEN
3457 nelem(l, ii) = nn
3458 EXIT
3459 ELSE
3460 nelem(l, ii) = ll
3461 nn = nn - ll
3462 END IF
3463 END DO
3464 END DO
3465 ncalc = nelem - ncore
3466 END IF
3467
3468 ! readjust the occupation number of the orbitals as requested by user
3469 ! this is done to break symmetry (bs) and bias the initial guess
3470 ! to the pre-defined multiplicity/charge state of the atom
3471 IF (bs_occupation) THEN
3472 DO is = 1, 2
3473 DO i = 1, SIZE(addel, 1)
3474 ne = addel(i, is)
3475 l = laddel(i, is)
3476 nn = naddel(i, is) - l
3477 IF (ne /= 0) THEN
3478 IF (nn == 0) THEN
3479 DO ii = SIZE(nelem, 2), 1, -1
3480 IF (ncalc(l, ii) > 0) THEN
3481 IF ((ncalc(l, ii) + ne) < 2*(2*l + 1) + 1) THEN
3482 edelta(l, ii, is) = edelta(l, ii, is) + ne
3483 nn = ii
3484 ELSE
3485 edelta(l, ii + 1, is) = edelta(l, ii + 1, is) + ne
3486 nn = ii + 1
3487 END IF
3488 EXIT
3489 ELSE IF (ii == 1) THEN
3490 edelta(l, ii, is) = edelta(l, ii, is) + ne
3491 nn = ii
3492 END IF
3493 END DO
3494 ELSE
3495 edelta(l, nn, is) = edelta(l, nn, is) + ne
3496 END IF
3497 IF (ncalc(l, nn) + edelta(l, nn, is) < 0) THEN
3498 edelta(l, nn, is) = -ncalc(l, nn)
3499 END IF
3500 END IF
3501 END DO
3502 END DO
3503 edelta = 0.5_dp*edelta
3504 ELSE IF (magnetization /= 0.0_dp) THEN
3505 dmag = 0.5_dp*abs(magnetization)
3506 DO l = 0, min(lmat, ubound(ptable(z)%e_conv, 1))
3507 ll = 2*(2*l + 1)
3508 ii = 0
3509 DO i = 1, SIZE(ncalc, 2)
3510 IF (ncalc(l, i) == 0) cycle
3511 IF (ncalc(l, i) == ll) cycle
3512 IF (ncalc(l, i) > dmag .AND. (ll - ncalc(l, i)) > dmag) THEN
3513 ii = i
3514 EXIT
3515 END IF
3516 END DO
3517 IF (ii /= 0) THEN
3518 edelta(l, ii, 1) = magnetization*0.5_dp
3519 edelta(l, ii, 2) = -magnetization*0.5_dp
3520 EXIT
3521 END IF
3522 END DO
3523 IF (ii == 0) THEN
3524 CALL cp_abort(__location__, &
3525 "Magnetization value cannot be imposed for this atom type")
3526 END IF
3527 END IF
3528
3529 IF (qs_kind%ghost .OR. qs_kind%floating) THEN
3530 nelem = 0
3531 ncore = 0
3532 ncalc = 0
3533 edelta = 0.0_dp
3534 END IF
3535
3536 END SUBROUTINE init_atom_electronic_state
3537
3538! **************************************************************************************************
3539!> \brief ...
3540!> \param econf ...
3541!> \param z ...
3542!> \param ncalc ...
3543!> \param ncore ...
3544!> \param nelem ...
3545! **************************************************************************************************
3546 SUBROUTINE set_pseudo_state(econf, z, ncalc, ncore, nelem)
3547 INTEGER, DIMENSION(:), POINTER :: econf
3548 INTEGER, INTENT(IN) :: z
3549 INTEGER, DIMENSION(0:lmat, 10), INTENT(OUT) :: ncalc, ncore, nelem
3550
3551 CHARACTER(LEN=default_string_length) :: message
3552 INTEGER :: ii, iounit, l, ll, lmin, nc, nn
3553 INTEGER, DIMENSION(0:lmat) :: econfx
3554 TYPE(cp_logger_type), POINTER :: logger
3555
3556 NULLIFY (logger)
3557 logger => cp_get_default_logger()
3558 iounit = cp_logger_get_default_io_unit(logger)
3559
3560 econfx = 0
3561 econfx(0:SIZE(econf) - 1) = econf
3562 IF (sum(econf) >= 0) THEN
3563 lmin = min(lmat, ubound(ptable(z)%e_conv, 1))
3564 ! number of core electrons
3565 nc = z - sum(econf)
3566 ! setup ncore
3567 ncore = 0
3568 SELECT CASE (nc)
3569 CASE (0)
3570 CASE (2)
3571 ncore(0, 1) = 2
3572 CASE (10)
3573 ncore(0, 1) = 2
3574 ncore(0, 2) = 2
3575 ncore(1, 1) = 6
3576 CASE (18)
3577 ncore(0, 1) = 2
3578 ncore(0, 2) = 2
3579 ncore(0, 3) = 2
3580 ncore(1, 1) = 6
3581 ncore(1, 2) = 6
3582 CASE (28)
3583 ncore(0, 1) = 2
3584 ncore(0, 2) = 2
3585 ncore(0, 3) = 2
3586 ncore(1, 1) = 6
3587 ncore(1, 2) = 6
3588 ncore(2, 1) = 10
3589 CASE (36)
3590 ncore(0, 1) = 2
3591 ncore(0, 2) = 2
3592 ncore(0, 3) = 2
3593 ncore(0, 4) = 2
3594 ncore(1, 1) = 6
3595 ncore(1, 2) = 6
3596 ncore(1, 3) = 6
3597 ncore(2, 1) = 10
3598 CASE (46)
3599 ncore(0, 1) = 2
3600 ncore(0, 2) = 2
3601 ncore(0, 3) = 2
3602 ncore(0, 4) = 2
3603 ncore(1, 1) = 6
3604 ncore(1, 2) = 6
3605 ncore(1, 3) = 6
3606 ncore(2, 1) = 10
3607 ncore(2, 2) = 10
3608 CASE (54)
3609 ncore(0, 1) = 2
3610 ncore(0, 2) = 2
3611 ncore(0, 3) = 2
3612 ncore(0, 4) = 2
3613 ncore(0, 5) = 2
3614 ncore(1, 1) = 6
3615 ncore(1, 2) = 6
3616 ncore(1, 3) = 6
3617 ncore(1, 4) = 6
3618 ncore(2, 1) = 10
3619 ncore(2, 2) = 10
3620 CASE (60)
3621 ncore(0, 1) = 2
3622 ncore(0, 2) = 2
3623 ncore(0, 3) = 2
3624 ncore(0, 4) = 2
3625 ncore(1, 1) = 6
3626 ncore(1, 2) = 6
3627 ncore(1, 3) = 6
3628 ncore(2, 1) = 10
3629 ncore(2, 2) = 10
3630 ncore(3, 1) = 14
3631 CASE (68)
3632 ncore(0, 1) = 2
3633 ncore(0, 2) = 2
3634 ncore(0, 3) = 2
3635 ncore(0, 4) = 2
3636 ncore(0, 5) = 2
3637 ncore(1, 1) = 6
3638 ncore(1, 2) = 6
3639 ncore(1, 3) = 6
3640 ncore(1, 4) = 6
3641 ncore(2, 1) = 10
3642 ncore(2, 2) = 10
3643 ncore(3, 1) = 14
3644 CASE (78)
3645 ncore(0, 1) = 2
3646 ncore(0, 2) = 2
3647 ncore(0, 3) = 2
3648 ncore(0, 4) = 2
3649 ncore(0, 5) = 2
3650 ncore(1, 1) = 6
3651 ncore(1, 2) = 6
3652 ncore(1, 3) = 6
3653 ncore(1, 4) = 6
3654 ncore(2, 1) = 10
3655 ncore(2, 2) = 10
3656 ncore(2, 3) = 10
3657 ncore(3, 1) = 14
3658 ! 79 - 92 5f incore PP
3659 CASE (79)
3660 ncore(0, 1) = 2
3661 ncore(0, 2) = 2
3662 ncore(0, 3) = 2
3663 ncore(0, 4) = 2
3664 ncore(0, 5) = 2
3665 ncore(1, 1) = 6
3666 ncore(1, 2) = 6
3667 ncore(1, 3) = 6
3668 ncore(1, 4) = 6
3669 ncore(2, 1) = 10
3670 ncore(2, 2) = 10
3671 ncore(2, 3) = 10
3672 ncore(3, 1) = 14
3673 ncore(3, 2) = 1
3674 CASE (80)
3675 ncore(0, 1) = 2
3676 ncore(0, 2) = 2
3677 ncore(0, 3) = 2
3678 ncore(0, 4) = 2
3679 ncore(0, 5) = 2
3680 ncore(1, 1) = 6
3681 ncore(1, 2) = 6
3682 ncore(1, 3) = 6
3683 ncore(1, 4) = 6
3684 ncore(2, 1) = 10
3685 ncore(2, 2) = 10
3686 ncore(2, 3) = 10
3687 ncore(3, 1) = 14
3688 ncore(3, 2) = 2
3689 CASE (81)
3690 ncore(0, 1) = 2
3691 ncore(0, 2) = 2
3692 ncore(0, 3) = 2
3693 ncore(0, 4) = 2
3694 ncore(0, 5) = 2
3695 ncore(1, 1) = 6
3696 ncore(1, 2) = 6
3697 ncore(1, 3) = 6
3698 ncore(1, 4) = 6
3699 ncore(2, 1) = 10
3700 ncore(2, 2) = 10
3701 ncore(2, 3) = 10
3702 ncore(3, 1) = 14
3703 ncore(3, 2) = 3
3704 CASE (82)
3705 ncore(0, 1) = 2
3706 ncore(0, 2) = 2
3707 ncore(0, 3) = 2
3708 ncore(0, 4) = 2
3709 ncore(0, 5) = 2
3710 ncore(1, 1) = 6
3711 ncore(1, 2) = 6
3712 ncore(1, 3) = 6
3713 ncore(1, 4) = 6
3714 ncore(2, 1) = 10
3715 ncore(2, 2) = 10
3716 ncore(2, 3) = 10
3717 ncore(3, 1) = 14
3718 ncore(3, 2) = 4
3719 CASE (83)
3720 ncore(0, 1) = 2
3721 ncore(0, 2) = 2
3722 ncore(0, 3) = 2
3723 ncore(0, 4) = 2
3724 ncore(0, 5) = 2
3725 ncore(1, 1) = 6
3726 ncore(1, 2) = 6
3727 ncore(1, 3) = 6
3728 ncore(1, 4) = 6
3729 ncore(2, 1) = 10
3730 ncore(2, 2) = 10
3731 ncore(2, 3) = 10
3732 ncore(3, 1) = 14
3733 ncore(3, 2) = 5
3734 CASE (84)
3735 ncore(0, 1) = 2
3736 ncore(0, 2) = 2
3737 ncore(0, 3) = 2
3738 ncore(0, 4) = 2
3739 ncore(0, 5) = 2
3740 ncore(1, 1) = 6
3741 ncore(1, 2) = 6
3742 ncore(1, 3) = 6
3743 ncore(1, 4) = 6
3744 ncore(2, 1) = 10
3745 ncore(2, 2) = 10
3746 ncore(2, 3) = 10
3747 ncore(3, 1) = 14
3748 ncore(3, 2) = 6
3749 CASE (85)
3750 ncore(0, 1) = 2
3751 ncore(0, 2) = 2
3752 ncore(0, 3) = 2
3753 ncore(0, 4) = 2
3754 ncore(0, 5) = 2
3755 ncore(1, 1) = 6
3756 ncore(1, 2) = 6
3757 ncore(1, 3) = 6
3758 ncore(1, 4) = 6
3759 ncore(2, 1) = 10
3760 ncore(2, 2) = 10
3761 ncore(2, 3) = 10
3762 ncore(3, 1) = 14
3763 ncore(3, 2) = 7
3764 CASE (86)
3765 ! this is not Rn core, add double assignment below
3766 ncore(0, 1) = 2
3767 ncore(0, 2) = 2
3768 ncore(0, 3) = 2
3769 ncore(0, 4) = 2
3770 ncore(0, 5) = 2
3771 ncore(1, 1) = 6
3772 ncore(1, 2) = 6
3773 ncore(1, 3) = 6
3774 ncore(1, 4) = 6
3775 ncore(2, 1) = 10
3776 ncore(2, 2) = 10
3777 ncore(2, 3) = 10
3778 ncore(3, 1) = 14
3779 ncore(3, 2) = 8
3780 CASE (87)
3781 ncore(0, 1) = 2
3782 ncore(0, 2) = 2
3783 ncore(0, 3) = 2
3784 ncore(0, 4) = 2
3785 ncore(0, 5) = 2
3786 ncore(1, 1) = 6
3787 ncore(1, 2) = 6
3788 ncore(1, 3) = 6
3789 ncore(1, 4) = 6
3790 ncore(2, 1) = 10
3791 ncore(2, 2) = 10
3792 ncore(2, 3) = 10
3793 ncore(3, 1) = 14
3794 ncore(3, 2) = 9
3795 CASE (88)
3796 ncore(0, 1) = 2
3797 ncore(0, 2) = 2
3798 ncore(0, 3) = 2
3799 ncore(0, 4) = 2
3800 ncore(0, 5) = 2
3801 ncore(1, 1) = 6
3802 ncore(1, 2) = 6
3803 ncore(1, 3) = 6
3804 ncore(1, 4) = 6
3805 ncore(2, 1) = 10
3806 ncore(2, 2) = 10
3807 ncore(2, 3) = 10
3808 ncore(3, 1) = 14
3809 ncore(3, 2) = 10
3810 CASE (89)
3811 ncore(0, 1) = 2
3812 ncore(0, 2) = 2
3813 ncore(0, 3) = 2
3814 ncore(0, 4) = 2
3815 ncore(0, 5) = 2
3816 ncore(1, 1) = 6
3817 ncore(1, 2) = 6
3818 ncore(1, 3) = 6
3819 ncore(1, 4) = 6
3820 ncore(2, 1) = 10
3821 ncore(2, 2) = 10
3822 ncore(2, 3) = 10
3823 ncore(3, 1) = 14
3824 ncore(3, 2) = 11
3825 CASE (90)
3826 ncore(0, 1) = 2
3827 ncore(0, 2) = 2
3828 ncore(0, 3) = 2
3829 ncore(0, 4) = 2
3830 ncore(0, 5) = 2
3831 ncore(1, 1) = 6
3832 ncore(1, 2) = 6
3833 ncore(1, 3) = 6
3834 ncore(1, 4) = 6
3835 ncore(2, 1) = 10
3836 ncore(2, 2) = 10
3837 ncore(2, 3) = 10
3838 ncore(3, 1) = 14
3839 ncore(3, 2) = 12
3840 CASE (91)
3841 ncore(0, 1) = 2
3842 ncore(0, 2) = 2
3843 ncore(0, 3) = 2
3844 ncore(0, 4) = 2
3845 ncore(0, 5) = 2
3846 ncore(1, 1) = 6
3847 ncore(1, 2) = 6
3848 ncore(1, 3) = 6
3849 ncore(1, 4) = 6
3850 ncore(2, 1) = 10
3851 ncore(2, 2) = 10
3852 ncore(2, 3) = 10
3853 ncore(3, 1) = 14
3854 ncore(3, 2) = 13
3855 CASE (92)
3856 ncore(0, 1) = 2
3857 ncore(0, 2) = 2
3858 ncore(0, 3) = 2
3859 ncore(0, 4) = 2
3860 ncore(0, 5) = 2
3861 ncore(1, 1) = 6
3862 ncore(1, 2) = 6
3863 ncore(1, 3) = 6
3864 ncore(1, 4) = 6
3865 ncore(2, 1) = 10
3866 ncore(2, 2) = 10
3867 ncore(2, 3) = 10
3868 ncore(3, 1) = 14
3869 ncore(3, 2) = 14
3870 CASE DEFAULT
3871 ncore(0, 1) = -1
3872 END SELECT
3873 ! special cases of double assignments
3874 IF (z == 65 .AND. econfx(3) == 0) THEN
3875 ! 4f in core for Tb
3876 ncore = 0
3877 ncore(0, 1) = -1
3878 END IF
3879 ! if there is still no core, check for special cases
3880 IF (ncore(0, 1) <= 0) THEN
3881 IF (z >= 58 .AND. z <= 71) THEN
3882 ! 4f-in-core PPs for lanthanides
3883 nc = z - sum(econf)
3884 ! setup ncore
3885 ncore = 0
3886 SELECT CASE (nc)
3887 CASE (29:42)
3888 ncore(0, 1) = 2
3889 ncore(0, 2) = 2
3890 ncore(0, 3) = 2
3891 ncore(1, 1) = 6
3892 ncore(1, 2) = 6
3893 ncore(2, 1) = 10
3894 ncore(3, 1) = nc - 28
3895 message = "A small-core pseudopotential with 4f-in-core is used for the lanthanide "// &
3896 trim(ptable(z)%symbol)
3897 cphint(trim(message))
3898 CASE (47:60)
3899 ncore(0, 1) = 2
3900 ncore(0, 2) = 2
3901 ncore(0, 3) = 2
3902 ncore(0, 4) = 2
3903 ncore(1, 1) = 6
3904 ncore(1, 2) = 6
3905 ncore(1, 3) = 6
3906 ncore(2, 1) = 10
3907 ncore(2, 2) = 10
3908 ncore(3, 1) = nc - 46
3909 message = "A medium-core pseudopotential with 4f-in-core is used for the lanthanide "// &
3910 trim(ptable(z)%symbol)
3911 cphint(trim(message))
3912 CASE DEFAULT
3913 ncore(0, 1) = -1
3914 END SELECT
3915 END IF
3916 END IF
3917 ! if the core is established, finish the setup
3918 IF (ncore(0, 1) >= 0) THEN
3919 DO l = 0, lmin
3920 ll = 2*(2*l + 1)
3921 nn = sum(ncore(l, :)) + econfx(l)
3922 ii = 0
3923 DO
3924 ii = ii + 1
3925 IF (nn <= ll) THEN
3926 nelem(l, ii) = nn
3927 EXIT
3928 ELSE
3929 nelem(l, ii) = ll
3930 nn = nn - ll
3931 END IF
3932 END DO
3933 END DO
3934 ncalc = nelem - ncore
3935 ELSE
3936 ! test for compatibility of valence occupation and full atomic occupation
3937 IF (iounit > 0) THEN
3938 WRITE (iounit, "(/,A,A2)") "WARNING: Core states irregular for atom type ", ptable(z)%symbol
3939 WRITE (iounit, "(A,10I3)") "WARNING: Redefine ELEC_CONF in the KIND section"
3940 cpabort("Incompatible Atomic Occupations Detected")
3941 END IF
3942 END IF
3943 ELSE
3944 lmin = min(lmat, ubound(ptable(z)%e_conv, 1))
3945 ncore = 0
3946 ncalc = 0
3947 DO l = 0, lmin
3948 ll = 2*(2*l + 1)
3949 nn = abs(econfx(l))
3950 ii = 0
3951 DO
3952 ii = ii + 1
3953 IF (nn <= ll) THEN
3954 ncalc(l, ii) = -nn
3955 EXIT
3956 ELSE
3957 ncalc(l, ii) = -ll
3958 nn = nn - ll
3959 END IF
3960 END DO
3961 END DO
3962 nelem = ncalc
3963 END IF
3964
3965 END SUBROUTINE set_pseudo_state
3966
3967! **************************************************************************************************
3968!> \brief finds if a given qs run needs to use nlcc
3969!> \param qs_kind_set ...
3970!> \return ...
3971! **************************************************************************************************
3972 FUNCTION has_nlcc(qs_kind_set) RESULT(nlcc)
3973
3974 TYPE(qs_kind_type), DIMENSION(:) :: qs_kind_set
3975 LOGICAL :: nlcc
3976
3977 INTEGER :: ikind
3978 LOGICAL :: nlcc_present
3979 TYPE(gth_potential_type), POINTER :: gth_potential
3980 TYPE(sgp_potential_type), POINTER :: sgp_potential
3981
3982 nlcc = .false.
3983
3984 DO ikind = 1, SIZE(qs_kind_set)
3985 CALL get_qs_kind(qs_kind_set(ikind), gth_potential=gth_potential, sgp_potential=sgp_potential)
3986 IF (ASSOCIATED(gth_potential)) THEN
3987 CALL get_potential(potential=gth_potential, nlcc_present=nlcc_present)
3988 nlcc = nlcc .OR. nlcc_present
3989 ELSEIF (ASSOCIATED(sgp_potential)) THEN
3990 CALL get_potential(potential=sgp_potential, has_nlcc=nlcc_present)
3991 nlcc = nlcc .OR. nlcc_present
3992 END IF
3993 END DO
3994
3995 END FUNCTION has_nlcc
3996
3997! **************************************************************************************************
3998
3999END MODULE qs_kind_types
static int imax(int x, int y)
Returns the larger of two given integers (missing from the C standard)
subroutine, public atom_sgp_release(sgp_pot)
...
Definition atom_sgp.F:767
subroutine, public sgp_construction(sgp_pot, ecp_pot, upf_pot, orb_basis, error)
...
Definition atom_sgp.F:79
Define the atom type and its sub types.
Definition atom_types.F:15
integer, parameter, public lmat
Definition atom_types.F:67
subroutine, public read_ecp_potential(element_symbol, potential, pseudo_name, pseudo_file, potential_section)
...
Routines that process Quantum Espresso UPF files.
Definition atom_upf.F:14
subroutine, public atom_read_upf(pot, upf_filename, read_header)
...
Definition atom_upf.F:102
pure subroutine, public atom_release_upf(upfpot)
...
Definition atom_upf.F:875
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 remove_basis_set_container(basis)
...
subroutine, public remove_basis_from_container(container, inum, basis_type)
...
subroutine, public get_basis_from_container(container, basis_set, inumbas, basis_type)
Retrieve a basis set from the container.
subroutine, public add_basis_set_to_container(container, basis_set, basis_set_type)
...
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)
...
subroutine, public write_orb_basis_set(orb_basis_set, output_unit, header)
Write a Gaussian-type orbital (GTO) basis set data set to the output unit.
subroutine, public deallocate_sto_basis_set(sto_basis_set)
...
subroutine, public init_aux_basis_set(gto_basis_set)
...
subroutine, public allocate_gto_basis_set(gto_basis_set)
...
subroutine, public combine_basis_sets(basis_set, basis_set_add)
...
subroutine, public write_gto_basis_set(gto_basis_set, output_unit, header)
Write a Gaussian-type orbital (GTO) basis set data set to the output unit.
subroutine, public allocate_sto_basis_set(sto_basis_set)
...
subroutine, public create_gto_from_sto_basis(sto_basis_set, gto_basis_set, ngauss, ortho)
...
subroutine, public read_sto_basis_set(element_symbol, basis_set_name, sto_basis_set, para_env, dft_section)
...
subroutine, public init_orb_basis_set(gto_basis_set)
Initialise a 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 ...
integer function, public cp_logger_get_default_io_unit(logger)
returns the unit nr for the ionode (-1 on all other processors) skips as well checks if the procs cal...
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...
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)
...
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, parameter, public cp_p_file
integer function, public cp_print_key_should_output(iteration_info, basis_section, print_key_path, used_print_key, first_time)
returns what should be done with the given property if btest(res,cp_p_store) then the property should...
Add the DFT+U contribution to the Hamiltonian matrix.
Definition dft_plus_u.F:18
Definition of the atomic potential types.
subroutine, public set_default_all_potential(potential, z, zeff_correction)
...
subroutine, public create_1c_basis(orb_basis, soft_basis, gapw_1c_basis, basis_1c_level)
create the one center basis from the orbital basis
collects all constants needed in input so that they can be used without circular dependencies
integer, parameter, public gapw_1c_large
integer, parameter, public do_method_pdg
integer, parameter, public do_method_pnnl
integer, parameter, public gapw_1c_medium
integer, parameter, public do_method_pw
integer, parameter, public do_method_rm1
integer, parameter, public gapw_1c_small
integer, parameter, public do_method_pm3
integer, parameter, public do_sirius
integer, parameter, public do_method_mndo
integer, parameter, public gapw_1c_orb
integer, parameter, public gapw_1c_very_large
integer, parameter, public do_method_mndod
integer, parameter, public do_method_am1
integer, parameter, public do_method_dftb
integer, parameter, public do_qs
integer, parameter, public do_method_xtb
integer, parameter, public do_method_pm6fm
integer, parameter, public do_method_pm6
objects that represent the structure of input sections and the data contained in an input section
recursive type(section_vals_type) function, pointer, public section_vals_get_subs_vals(section_vals, subsection_name, i_rep_section, can_return_null)
returns the values of the requested subsection
subroutine, public section_vals_get(section_vals, ref_count, n_repetition, n_subs_vals_rep, section, explicit)
returns various attributes about the section_vals
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
integer, parameter, public default_path_length
Definition kinds.F:58
Definition of mathematical constants and functions.
real(kind=dp), parameter, public pi
Interface to the message passing library MPI.
Provides Cartesian and spherical orbital pointers and indices.
subroutine, public init_orbital_pointers(maxl)
Initialize or update the orbital pointers.
integer, dimension(:), allocatable, public nco
integer, dimension(:), allocatable, public ncoset
Factory routines for potentials used e.g. by pao_param_exp and pao_ml.
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 allocate_paw_proj_set(paw_proj_set)
Allocate projector type for GAPW.
subroutine, public projectors(paw_proj, basis_1c, orb_basis, rc, qs_control, max_rad_local_type, force_env_section)
Initialize the projector-type set data set.
subroutine, public deallocate_paw_proj_set(paw_proj_set)
Deallocate a projector-type set data set.
Periodic Table related data definitions.
type(atom), dimension(0:nelem), public ptable
integer, parameter, public nelem
subroutine, public get_ptable_info(symbol, number, amass, ielement, covalent_radius, metallic_radius, vdw_radius, found)
Pass information about the kind given the element symbol.
Definition of physical constants:
Definition physcon.F:68
real(kind=dp), parameter, public evolt
Definition physcon.F:183
real(kind=dp), parameter, public angstrom
Definition physcon.F:144
real(kind=dp), parameter, public bohr
Definition physcon.F:147
Types used by CNEO-DFT (see J. Chem. Theory Comput. 2025, 21, 16, 7865–7877)
subroutine, public allocate_cneo_potential(potential)
...
subroutine, public write_cneo_potential(potential, output_unit)
...
subroutine, public deallocate_cneo_potential(potential)
...
subroutine, public set_cneo_potential(potential, z, mass, elec_conf, nsgf, nne, npsgf, nsotot, my_gcc_h, my_gcc_s, ovlp, kin, utrans, distance, harmonics, qlm_gg, gg, vgg, n2oindex, o2nindex, rad2l, oorad2l)
...
subroutine, public get_cneo_potential(potential, z, zeff, mass, elec_conf, nsgf, nne, npsgf, nsotot, my_gcc_h, my_gcc_s, ovlp, kin, utrans, distance, harmonics, qlm_gg, gg, vgg, n2oindex, o2nindex, rad2l, oorad2l)
...
Definition of the DFTB parameter types.
Working with the DFTB parameter types.
subroutine, public deallocate_dftb_atom_param(dftb_parameter)
...
subroutine, public write_dftb_atom_param(dftb_parameter, subsys_section)
...
subroutine, public get_dftb_atom_param(dftb_parameter, name, typ, defined, z, zeff, natorb, lmax, skself, occupation, eta, energy, cutoff, xi, di, rcdisp, dudq)
...
Definition of disperson types for DFT calculations.
subroutine, public deallocate_grid_atom(grid_atom)
Deallocate a Gaussian-type orbital (GTO) basis set data set.
subroutine, public allocate_grid_atom(grid_atom)
Initialize components of the grid_atom_type structure.
subroutine, public allocate_harmonics_atom(harmonics)
Allocate a spherical harmonics set for the atom grid.
subroutine, public deallocate_harmonics_atom(harmonics)
Deallocate the spherical harmonics set for the atom grid.
Define the quickstep kind type and their sub types.
subroutine, public set_pseudo_state(econf, z, ncalc, ncore, nelem)
...
logical function, public has_nlcc(qs_kind_set)
finds if a given qs run needs to use nlcc
subroutine, public deallocate_qs_kind_set(qs_kind_set)
Destructor routine for a set of qs kinds.
subroutine, public get_qs_kind(qs_kind, basis_set, basis_type, ncgf, nsgf, all_potential, tnadd_potential, gth_potential, sgp_potential, upf_potential, cneo_potential, se_parameter, dftb_parameter, xtb_parameter, dftb3_param, zatom, zeff, elec_conf, mao, lmax_dftb, alpha_core_charge, ccore_charge, core_charge, core_charge_radius, paw_proj_set, paw_atom, hard_radius, hard0_radius, max_rad_local, covalent_radius, vdw_radius, gpw_type_forced, harmonics, max_iso_not0, max_s_harm, grid_atom, ngrid_ang, ngrid_rad, lmax_rho0, dft_plus_u_atom, l_of_dft_plus_u, n_of_dft_plus_u, u_minus_j, u_of_dft_plus_u, j_of_dft_plus_u, alpha_of_dft_plus_u, beta_of_dft_plus_u, j0_of_dft_plus_u, occupation_of_dft_plus_u, dispersion, bs_occupation, magnetization, no_optimize, addel, laddel, naddel, orbitals, max_scf, eps_scf, smear, u_ramping, u_minus_j_target, eps_u_ramping, init_u_ramping_each_scf, reltmat, ghost, monovalent, floating, name, element_symbol, pao_basis_size, pao_model_file, pao_potentials, pao_descriptors, nelec)
Get attributes of an atomic kind.
subroutine, public init_atom_electronic_state(atomic_kind, qs_kind, ncalc, ncore, nelem, edelta)
...
subroutine, public get_qs_kind_set(qs_kind_set, all_potential_present, tnadd_potential_present, gth_potential_present, sgp_potential_present, paw_atom_present, dft_plus_u_atom_present, maxcgf, maxsgf, maxco, maxco_proj, maxgtops, maxlgto, maxlprj, maxnset, maxsgf_set, ncgf, npgf, nset, nsgf, nshell, maxpol, maxlppl, maxlppnl, maxppnl, nelectron, maxder, max_ngrid_rad, max_sph_harm, maxg_iso_not0, lmax_rho0, basis_rcut, basis_type, total_zeff_corr, npgf_seg, cneo_potential_present, nkind_q, natom_q)
Get attributes of an atomic kind set.
subroutine, public init_gapw_nlcc(qs_kind_set)
...
subroutine, public init_cneo_basis_set(qs_kind_set, qs_control)
...
subroutine, public set_qs_kind(qs_kind, paw_atom, ghost, floating, hard_radius, hard0_radius, covalent_radius, vdw_radius, lmax_rho0, zeff, no_optimize, dispersion, u_minus_j, reltmat, dftb_parameter, xtb_parameter, elec_conf, pao_basis_size)
Set the components of an atomic kind data set.
subroutine, public write_qs_kind_set(qs_kind_set, subsys_section)
Write an atomic kind set data set to the output unit.
subroutine, public init_gapw_basis_set(qs_kind_set, qs_control, force_env_section, modify_qs_control)
...
subroutine, public init_qs_kind_set(qs_kind_set)
Initialise an atomic kind set data set.
subroutine, public check_qs_kind_set(qs_kind_set, dft_control, subsys_section)
...
subroutine, public write_gto_basis_sets(qs_kind_set, subsys_section)
Write all the GTO basis sets of an atomic kind set to the output unit (for the printing of the unnorm...
subroutine, public create_qs_kind_set(qs_kind_set, atomic_kind_set, kind_section, para_env, force_env_section, silent)
Read an atomic kind set data set from the input file.
Definition of the semi empirical parameter types.
subroutine, public write_se_param(sep, subsys_section)
Writes the semi-empirical type.
subroutine, public semi_empirical_create(sep)
Allocate semi-empirical type.
subroutine, public get_se_param(sep, name, typ, defined, z, zeff, natorb, eheat, beta, sto_exponents, uss, upp, udd, uff, alp, eisol, gss, gsp, gpp, gp2, acoul, nr, de, ass, asp, app, hsp, gsd, gpd, gdd, ppddg, dpddg, ngauss)
Get info from the semi-empirical type.
subroutine, public semi_empirical_release(sep)
Deallocate the semi-empirical type.
Working with the semi empirical parameter types.
subroutine, public se_param_set_default(sep, z, method)
Initialize parameter for a semi_empirival type.
subroutine, public init_se_param(sep, orb_basis_set, ngauss)
Initialize semi_empirical type.
subroutine, public create_soft_basis(orb_basis, soft_basis, eps_fit, rc, paw_atom, paw_type_forced, gpw_r3d_rs_type_forced)
create the soft basis from a GTO basis
Utilities for string manipulations.
elemental subroutine, public uppercase(string)
Convert all lower case characters in a string to upper case.
Read xTB parameters.
real(kind=dp) function, public xtb_set_kab(za, zb, xtb_control)
...
Definition of the xTB parameter types.
Definition xtb_types.F:20
subroutine, public get_xtb_atom_param(xtb_parameter, symbol, aname, typ, defined, z, zeff, natorb, lmax, nao, lao, rcut, rcov, kx, eta, xgamma, alpha, zneff, nshell, nval, lval, kpoly, kappa, hen, zeta, xi, kappa0, alpg, occupation, electronegativity, chmax, en, kqat2, kcn, kq)
...
Definition xtb_types.F:199
subroutine, public deallocate_xtb_atom_param(xtb_parameter)
...
Definition xtb_types.F:149
subroutine, public write_xtb_atom_param(xtb_parameter, gfn_type, subsys_section)
...
Definition xtb_types.F:364
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...
stores all the informations relevant to an mpi environment
Holds information about a PAO descriptor.
Holds information about a PAO potential.
Provides all information about a quickstep kind.