Line data Source code
1 : !--------------------------------------------------------------------------------------------------!
2 : ! CP2K: A general program to perform molecular dynamics simulations !
3 : ! Copyright 2000-2026 CP2K developers group <https://cp2k.org> !
4 : ! !
5 : ! SPDX-License-Identifier: GPL-2.0-or-later !
6 : !--------------------------------------------------------------------------------------------------!
7 :
8 : ! **************************************************************************************************
9 : !> \brief 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 : ! **************************************************************************************************
23 : MODULE qs_kind_types
24 : USE atom_sgp, ONLY: atom_sgp_potential_type,&
25 : atom_sgp_release,&
26 : sgp_construction
27 : USE atom_types, ONLY: atom_ecppot_type,&
28 : lmat,&
29 : read_ecp_potential
30 : USE atom_upf, ONLY: atom_read_upf,&
31 : atom_release_upf,&
32 : atom_upfpot_type
33 : USE atomic_kind_types, ONLY: atomic_kind_type,&
34 : get_atomic_kind
35 : USE basis_set_container_types, ONLY: add_basis_set_to_container,&
36 : basis_set_container_type,&
37 : get_basis_from_container,&
38 : remove_basis_from_container,&
39 : remove_basis_set_container
40 : USE basis_set_types, ONLY: &
41 : allocate_gto_basis_set, allocate_sto_basis_set, combine_basis_sets, &
42 : create_gto_from_sto_basis, deallocate_sto_basis_set, get_gto_basis_set, &
43 : gto_basis_set_type, init_aux_basis_set, init_orb_basis_set, read_gto_basis_set, &
44 : read_sto_basis_set, sto_basis_set_type, write_gto_basis_set, write_orb_basis_set
45 : USE cp_control_types, ONLY: dft_control_type,&
46 : qs_control_type,&
47 : xtb_control_type
48 : USE cp_log_handling, ONLY: cp_get_default_logger,&
49 : cp_logger_get_default_io_unit,&
50 : cp_logger_type
51 : USE cp_output_handling, ONLY: cp_p_file,&
52 : cp_print_key_finished_output,&
53 : cp_print_key_should_output,&
54 : cp_print_key_unit_nr
55 : USE external_potential_types, ONLY: &
56 : all_potential_type, allocate_potential, deallocate_potential, get_potential, &
57 : gth_potential_type, init_potential, local_potential_type, read_potential, &
58 : set_default_all_potential, set_potential, sgp_potential_type, write_potential
59 : USE gapw_1c_basis_set, ONLY: create_1c_basis
60 : USE input_constants, ONLY: &
61 : do_method_am1, do_method_dftb, do_method_mndo, do_method_mndod, do_method_pdg, &
62 : do_method_pm3, do_method_pm6, do_method_pm6fm, do_method_pnnl, do_method_pw, &
63 : do_method_rm1, do_method_xtb, do_qs, do_sirius, gapw_1c_large, gapw_1c_medium, &
64 : gapw_1c_orb, gapw_1c_small, gapw_1c_very_large
65 : USE input_section_types, ONLY: section_vals_get,&
66 : section_vals_get_subs_vals,&
67 : section_vals_type,&
68 : section_vals_val_get
69 : USE kinds, ONLY: default_path_length,&
70 : default_string_length,&
71 : dp
72 : USE mathconstants, ONLY: pi
73 : USE message_passing, ONLY: mp_para_env_type
74 : USE orbital_pointers, ONLY: init_orbital_pointers,&
75 : nco,&
76 : ncoset
77 : USE paw_proj_set_types, ONLY: allocate_paw_proj_set,&
78 : deallocate_paw_proj_set,&
79 : get_paw_proj_set,&
80 : paw_proj_set_type,&
81 : projectors
82 : USE periodic_table, ONLY: get_ptable_info,&
83 : ptable
84 : USE physcon, ONLY: angstrom,&
85 : bohr,&
86 : evolt
87 : USE qs_cneo_types, ONLY: allocate_cneo_potential,&
88 : cneo_potential_type,&
89 : deallocate_cneo_potential,&
90 : get_cneo_potential,&
91 : set_cneo_potential,&
92 : write_cneo_potential
93 : USE qs_dftb_types, ONLY: qs_dftb_atom_type
94 : USE qs_dftb_utils, ONLY: deallocate_dftb_atom_param,&
95 : get_dftb_atom_param,&
96 : write_dftb_atom_param
97 : USE qs_dispersion_types, ONLY: qs_atom_dispersion_type
98 : USE qs_grid_atom, ONLY: allocate_grid_atom,&
99 : deallocate_grid_atom,&
100 : grid_atom_type
101 : USE qs_harmonics_atom, ONLY: allocate_harmonics_atom,&
102 : deallocate_harmonics_atom,&
103 : harmonics_atom_type
104 : USE semi_empirical_types, ONLY: get_se_param,&
105 : semi_empirical_create,&
106 : semi_empirical_release,&
107 : semi_empirical_type,&
108 : write_se_param
109 : USE semi_empirical_utils, ONLY: init_se_param,&
110 : se_param_set_default
111 : USE soft_basis_set, ONLY: create_soft_basis
112 : USE string_utilities, ONLY: uppercase
113 : USE xtb_parameters, ONLY: xtb_set_kab
114 : USE xtb_types, ONLY: deallocate_xtb_atom_param,&
115 : get_xtb_atom_param,&
116 : write_xtb_atom_param,&
117 : xtb_atom_type
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 : ! **************************************************************************************************
155 : TYPE pao_potential_type
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 : ! **************************************************************************************************
166 : TYPE pao_descriptor_type
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 : ! **************************************************************************************************
177 : TYPE qs_kind_type
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 : !
192 : TYPE(basis_set_container_type), &
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, &
246 : deallocate_qs_kind_set, &
247 : get_qs_kind, &
248 : get_qs_kind_set, &
249 : has_nlcc, &
250 : init_qs_kind_set, &
251 : init_gapw_basis_set, &
252 : init_gapw_nlcc, &
253 : create_qs_kind_set, &
254 : set_qs_kind, &
255 : write_qs_kind_set, &
256 : write_gto_basis_sets, &
257 : init_atom_electronic_state, set_pseudo_state, &
258 : init_cneo_basis_set
259 :
260 : ! Public data types
261 : PUBLIC :: qs_kind_type, pao_potential_type, pao_descriptor_type
262 :
263 : CONTAINS
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 8858 : 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 8858 : IF (ASSOCIATED(qs_kind_set)) THEN
279 :
280 8858 : nkind = SIZE(qs_kind_set)
281 :
282 25585 : DO ikind = 1, nkind
283 16727 : IF (ASSOCIATED(qs_kind_set(ikind)%all_potential)) THEN
284 6784 : CALL deallocate_potential(qs_kind_set(ikind)%all_potential)
285 : END IF
286 16727 : IF (ASSOCIATED(qs_kind_set(ikind)%tnadd_potential)) THEN
287 28 : CALL deallocate_potential(qs_kind_set(ikind)%tnadd_potential)
288 : END IF
289 16727 : IF (ASSOCIATED(qs_kind_set(ikind)%gth_potential)) THEN
290 9675 : CALL deallocate_potential(qs_kind_set(ikind)%gth_potential)
291 : END IF
292 16727 : IF (ASSOCIATED(qs_kind_set(ikind)%sgp_potential)) THEN
293 88 : CALL deallocate_potential(qs_kind_set(ikind)%sgp_potential)
294 : END IF
295 16727 : IF (ASSOCIATED(qs_kind_set(ikind)%upf_potential)) THEN
296 20 : CALL atom_release_upf(qs_kind_set(ikind)%upf_potential)
297 20 : DEALLOCATE (qs_kind_set(ikind)%upf_potential)
298 : END IF
299 16727 : IF (ASSOCIATED(qs_kind_set(ikind)%cneo_potential)) THEN
300 8 : CALL deallocate_cneo_potential(qs_kind_set(ikind)%cneo_potential)
301 : END IF
302 16727 : IF (ASSOCIATED(qs_kind_set(ikind)%se_parameter)) THEN
303 2244 : CALL semi_empirical_release(qs_kind_set(ikind)%se_parameter)
304 : END IF
305 16727 : IF (ASSOCIATED(qs_kind_set(ikind)%dftb_parameter)) THEN
306 618 : CALL deallocate_dftb_atom_param(qs_kind_set(ikind)%dftb_parameter)
307 : END IF
308 16727 : IF (ASSOCIATED(qs_kind_set(ikind)%xtb_parameter)) THEN
309 2562 : CALL deallocate_xtb_atom_param(qs_kind_set(ikind)%xtb_parameter)
310 : END IF
311 16727 : IF (ASSOCIATED(qs_kind_set(ikind)%paw_proj_set)) THEN
312 2376 : CALL deallocate_paw_proj_set(qs_kind_set(ikind)%paw_proj_set)
313 : END IF
314 16727 : IF (ASSOCIATED(qs_kind_set(ikind)%harmonics)) THEN
315 2788 : CALL deallocate_harmonics_atom(qs_kind_set(ikind)%harmonics)
316 : END IF
317 16727 : IF (ASSOCIATED(qs_kind_set(ikind)%grid_atom)) THEN
318 2788 : CALL deallocate_grid_atom(qs_kind_set(ikind)%grid_atom)
319 : END IF
320 16727 : IF (ASSOCIATED(qs_kind_set(ikind)%elec_conf)) THEN
321 16297 : DEALLOCATE (qs_kind_set(ikind)%elec_conf)
322 : END IF
323 :
324 16727 : IF (ASSOCIATED(qs_kind_set(ikind)%dft_plus_u)) THEN
325 36 : IF (ASSOCIATED(qs_kind_set(ikind)%dft_plus_u%orbitals)) THEN
326 4 : DEALLOCATE (qs_kind_set(ikind)%dft_plus_u%orbitals)
327 : END IF
328 36 : IF (ASSOCIATED(qs_kind_set(ikind)%dft_plus_u%nelec)) THEN
329 4 : DEALLOCATE (qs_kind_set(ikind)%dft_plus_u%nelec)
330 : END IF
331 36 : DEALLOCATE (qs_kind_set(ikind)%dft_plus_u)
332 : END IF
333 :
334 16727 : IF (ASSOCIATED(qs_kind_set(ikind)%nlcc_pot)) THEN
335 10 : DEALLOCATE (qs_kind_set(ikind)%nlcc_pot)
336 : END IF
337 :
338 16727 : IF (ASSOCIATED(qs_kind_set(ikind)%dispersion)) THEN
339 2728 : DEALLOCATE (qs_kind_set(ikind)%dispersion)
340 : END IF
341 16727 : IF (ASSOCIATED(qs_kind_set(ikind)%addel)) THEN
342 66 : DEALLOCATE (qs_kind_set(ikind)%addel)
343 : END IF
344 16727 : IF (ASSOCIATED(qs_kind_set(ikind)%naddel)) THEN
345 66 : DEALLOCATE (qs_kind_set(ikind)%naddel)
346 : END IF
347 16727 : IF (ASSOCIATED(qs_kind_set(ikind)%laddel)) THEN
348 66 : DEALLOCATE (qs_kind_set(ikind)%laddel)
349 : END IF
350 16727 : IF (ASSOCIATED(qs_kind_set(ikind)%reltmat)) THEN
351 28 : DEALLOCATE (qs_kind_set(ikind)%reltmat)
352 : END IF
353 :
354 16727 : IF (ASSOCIATED(qs_kind_set(ikind)%pao_potentials)) THEN
355 11237 : DEALLOCATE (qs_kind_set(ikind)%pao_potentials)
356 : END IF
357 16727 : IF (ASSOCIATED(qs_kind_set(ikind)%pao_descriptors)) THEN
358 11237 : DEALLOCATE (qs_kind_set(ikind)%pao_descriptors)
359 : END IF
360 :
361 25585 : CALL remove_basis_set_container(qs_kind_set(ikind)%basis_sets)
362 :
363 : END DO
364 8858 : DEALLOCATE (qs_kind_set)
365 : ELSE
366 : CALL cp_abort(__LOCATION__, &
367 : "The pointer qs_kind_set is not associated and "// &
368 0 : "cannot be deallocated")
369 : END IF
370 :
371 8858 : 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 77338383 : 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 77338383 : IF (PRESENT(basis_type)) THEN
538 18191628 : my_basis_type = basis_type
539 : ELSE
540 59146755 : my_basis_type = "ORB"
541 : END IF
542 :
543 77338383 : IF (PRESENT(basis_set)) THEN
544 : CALL get_basis_from_container(qs_kind%basis_sets, basis_set=basis_set, &
545 19059633 : basis_type=my_basis_type)
546 : END IF
547 :
548 77338383 : IF (PRESENT(ncgf)) THEN
549 : CALL get_basis_from_container(qs_kind%basis_sets, basis_set=tmp_basis_set, &
550 960 : basis_type=my_basis_type)
551 960 : IF (ASSOCIATED(tmp_basis_set)) THEN
552 960 : CALL get_gto_basis_set(gto_basis_set=tmp_basis_set, ncgf=ncgf)
553 0 : ELSE IF (ASSOCIATED(qs_kind%dftb_parameter)) THEN
554 0 : l = qs_kind%dftb_parameter%lmax
555 0 : ncgf = ((l + 1)*(l + 2)*(l + 3))/6
556 : ELSE
557 0 : ncgf = 0
558 : END IF
559 : END IF
560 :
561 77338383 : IF (PRESENT(nsgf)) THEN
562 : CALL get_basis_from_container(qs_kind%basis_sets, basis_set=tmp_basis_set, &
563 300931 : basis_type=my_basis_type)
564 300931 : IF (ASSOCIATED(tmp_basis_set)) THEN
565 187915 : CALL get_gto_basis_set(gto_basis_set=tmp_basis_set, nsgf=nsgf)
566 113016 : ELSE IF (ASSOCIATED(qs_kind%dftb_parameter)) THEN
567 113012 : nsgf = qs_kind%dftb_parameter%natorb
568 : ELSE
569 4 : nsgf = 0
570 : END IF
571 : END IF
572 :
573 77338383 : IF (PRESENT(all_potential)) all_potential => qs_kind%all_potential
574 77338383 : IF (PRESENT(tnadd_potential)) tnadd_potential => qs_kind%tnadd_potential
575 77338383 : IF (PRESENT(gth_potential)) gth_potential => qs_kind%gth_potential
576 77338383 : IF (PRESENT(sgp_potential)) sgp_potential => qs_kind%sgp_potential
577 77338383 : IF (PRESENT(upf_potential)) upf_potential => qs_kind%upf_potential
578 77338383 : IF (PRESENT(cneo_potential)) cneo_potential => qs_kind%cneo_potential
579 77338383 : IF (PRESENT(se_parameter)) se_parameter => qs_kind%se_parameter
580 77338383 : IF (PRESENT(dftb_parameter)) dftb_parameter => qs_kind%dftb_parameter
581 77338383 : IF (PRESENT(xtb_parameter)) xtb_parameter => qs_kind%xtb_parameter
582 77338383 : IF (PRESENT(element_symbol)) element_symbol = qs_kind%element_symbol
583 77338383 : IF (PRESENT(name)) name = qs_kind%name
584 77338383 : IF (PRESENT(dftb3_param)) dftb3_param = qs_kind%dudq_dftb3
585 77338383 : IF (PRESENT(elec_conf)) elec_conf => qs_kind%elec_conf
586 77338383 : IF (PRESENT(alpha_core_charge)) THEN
587 241371 : IF (ASSOCIATED(qs_kind%all_potential)) THEN
588 : CALL get_potential(potential=qs_kind%all_potential, &
589 64894 : alpha_core_charge=alpha_core_charge)
590 176477 : ELSE IF (ASSOCIATED(qs_kind%gth_potential)) THEN
591 : CALL get_potential(potential=qs_kind%gth_potential, &
592 173927 : alpha_core_charge=alpha_core_charge)
593 2550 : ELSE IF (ASSOCIATED(qs_kind%sgp_potential)) THEN
594 : CALL get_potential(potential=qs_kind%sgp_potential, &
595 1134 : alpha_core_charge=alpha_core_charge)
596 1416 : ELSE IF (ASSOCIATED(qs_kind%cneo_potential)) THEN
597 0 : CPABORT("CNEO ALPHA CORE CHARGE NOT AVAILABLE")
598 : ELSE
599 1416 : alpha_core_charge = 1.0_dp
600 : END IF
601 : END IF
602 77338383 : IF (PRESENT(ccore_charge)) THEN
603 97375 : IF (ASSOCIATED(qs_kind%all_potential)) THEN
604 : CALL get_potential(potential=qs_kind%all_potential, &
605 15730 : ccore_charge=ccore_charge)
606 81645 : ELSE IF (ASSOCIATED(qs_kind%gth_potential)) THEN
607 : CALL get_potential(potential=qs_kind%gth_potential, &
608 80271 : ccore_charge=ccore_charge)
609 1374 : ELSE IF (ASSOCIATED(qs_kind%sgp_potential)) THEN
610 : CALL get_potential(potential=qs_kind%sgp_potential, &
611 560 : ccore_charge=ccore_charge)
612 814 : ELSE IF (ASSOCIATED(qs_kind%upf_potential)) THEN
613 0 : CPABORT("UPF CCORE CHARGE NOT AVAILABLE")
614 814 : ELSE IF (ASSOCIATED(qs_kind%cneo_potential)) THEN
615 0 : CPABORT("CNEO CCORE CHARGE NOT AVAILABLE")
616 : ELSE
617 814 : ccore_charge = 0.0_dp
618 : END IF
619 : END IF
620 77338383 : IF (PRESENT(core_charge_radius)) THEN
621 99845 : IF (ASSOCIATED(qs_kind%all_potential)) THEN
622 : CALL get_potential(potential=qs_kind%all_potential, &
623 44616 : core_charge_radius=core_charge_radius)
624 55229 : ELSE IF (ASSOCIATED(qs_kind%gth_potential)) THEN
625 : CALL get_potential(potential=qs_kind%gth_potential, &
626 54511 : core_charge_radius=core_charge_radius)
627 718 : ELSE IF (ASSOCIATED(qs_kind%sgp_potential)) THEN
628 : CALL get_potential(potential=qs_kind%sgp_potential, &
629 328 : core_charge_radius=core_charge_radius)
630 390 : ELSE IF (ASSOCIATED(qs_kind%upf_potential)) THEN
631 0 : CPABORT("UPF CORE CHARGE RADIUS NOT AVAILABLE")
632 390 : ELSE IF (ASSOCIATED(qs_kind%cneo_potential)) THEN
633 0 : CPABORT("CNEO CORE CHARGE RADIUS NOT AVAILABLE")
634 : ELSE
635 390 : core_charge_radius = 0.0_dp
636 : END IF
637 : END IF
638 77338383 : IF (PRESENT(core_charge)) THEN
639 47509 : IF (ASSOCIATED(qs_kind%all_potential)) THEN
640 : CALL get_potential(potential=qs_kind%all_potential, &
641 1064 : zeff=core_charge)
642 46445 : ELSE IF (ASSOCIATED(qs_kind%gth_potential)) THEN
643 : CALL get_potential(potential=qs_kind%gth_potential, &
644 46445 : zeff=core_charge)
645 0 : ELSE IF (ASSOCIATED(qs_kind%sgp_potential)) THEN
646 : CALL get_potential(potential=qs_kind%sgp_potential, &
647 0 : zeff=core_charge)
648 0 : ELSE IF (ASSOCIATED(qs_kind%upf_potential)) THEN
649 0 : CPABORT("UPF CORE CHARGE NOT AVAILABLE")
650 0 : ELSE IF (ASSOCIATED(qs_kind%cneo_potential)) THEN
651 : CALL get_cneo_potential(potential=qs_kind%cneo_potential, &
652 0 : zeff=core_charge)
653 : ELSE
654 0 : core_charge = 0.0_dp
655 : END IF
656 : END IF
657 :
658 77338383 : IF (PRESENT(zatom)) THEN
659 : ! Retrieve information on element
660 303692 : CALL get_ptable_info(qs_kind%element_symbol, ielement=zatom, found=found)
661 303692 : CPASSERT(found)
662 : END IF
663 :
664 77338383 : IF (PRESENT(zeff)) THEN
665 335318 : IF (ASSOCIATED(qs_kind%all_potential)) THEN
666 86645 : CALL get_potential(potential=qs_kind%all_potential, zeff=zeff)
667 248673 : ELSE IF (ASSOCIATED(qs_kind%gth_potential)) THEN
668 245697 : CALL get_potential(potential=qs_kind%gth_potential, zeff=zeff)
669 2976 : ELSE IF (ASSOCIATED(qs_kind%sgp_potential)) THEN
670 1611 : CALL get_potential(potential=qs_kind%sgp_potential, zeff=zeff)
671 1365 : ELSE IF (ASSOCIATED(qs_kind%upf_potential)) THEN
672 65 : zeff = qs_kind%upf_potential%zion
673 1300 : ELSE IF (ASSOCIATED(qs_kind%cneo_potential)) THEN
674 143 : CALL get_cneo_potential(potential=qs_kind%cneo_potential, zeff=zeff)
675 : ELSE
676 1157 : zeff = 0.0_dp
677 : END IF
678 : END IF
679 :
680 77338383 : IF (PRESENT(covalent_radius)) covalent_radius = qs_kind%covalent_radius
681 77338383 : IF (PRESENT(vdw_radius)) vdw_radius = qs_kind%vdw_radius
682 :
683 77338383 : IF (PRESENT(paw_proj_set)) paw_proj_set => qs_kind%paw_proj_set
684 77338383 : IF (PRESENT(paw_atom)) paw_atom = qs_kind%paw_atom
685 77338383 : IF (PRESENT(gpw_type_forced)) gpw_type_forced = qs_kind%gpw_type_forced
686 77338383 : IF (PRESENT(hard_radius)) hard_radius = qs_kind%hard_radius
687 77338383 : IF (PRESENT(hard0_radius)) hard0_radius = qs_kind%hard0_radius
688 77338383 : IF (PRESENT(max_rad_local)) max_rad_local = qs_kind%max_rad_local
689 77338383 : IF (PRESENT(harmonics)) harmonics => qs_kind%harmonics
690 77338383 : IF (PRESENT(max_s_harm)) THEN
691 9173114 : IF (ASSOCIATED(qs_kind%harmonics)) THEN
692 455018 : max_s_harm = qs_kind%harmonics%max_s_harm
693 : ELSE
694 8718096 : max_s_harm = 0
695 : END IF
696 : END IF
697 77338383 : IF (PRESENT(max_iso_not0)) THEN
698 9223614 : IF (ASSOCIATED(qs_kind%harmonics)) THEN
699 505518 : max_iso_not0 = qs_kind%harmonics%max_iso_not0
700 : ELSE
701 8718096 : max_iso_not0 = 0
702 : END IF
703 : END IF
704 77338383 : IF (PRESENT(grid_atom)) grid_atom => qs_kind%grid_atom
705 77338383 : IF (PRESENT(ngrid_ang)) ngrid_ang = qs_kind%ngrid_ang
706 77338383 : IF (PRESENT(ngrid_rad)) ngrid_rad = qs_kind%ngrid_rad
707 77338383 : IF (PRESENT(lmax_rho0)) lmax_rho0 = qs_kind%lmax_rho0
708 77338383 : IF (PRESENT(ghost)) ghost = qs_kind%ghost
709 77338383 : IF (PRESENT(monovalent)) monovalent = qs_kind%monovalent
710 77338383 : IF (PRESENT(floating)) floating = qs_kind%floating
711 77338383 : IF (PRESENT(dft_plus_u_atom)) dft_plus_u_atom = ASSOCIATED(qs_kind%dft_plus_u)
712 77338383 : IF (PRESENT(l_of_dft_plus_u)) THEN
713 5422 : IF (ASSOCIATED(qs_kind%dft_plus_u)) THEN
714 2698 : l_of_dft_plus_u = qs_kind%dft_plus_u%l
715 : ELSE
716 2724 : l_of_dft_plus_u = -1
717 : END IF
718 : END IF
719 77338383 : IF (PRESENT(n_of_dft_plus_u)) THEN
720 26 : IF (ASSOCIATED(qs_kind%dft_plus_u)) THEN
721 0 : n_of_dft_plus_u = qs_kind%dft_plus_u%n
722 : ELSE
723 26 : n_of_dft_plus_u = -1
724 : END IF
725 : END IF
726 77338383 : IF (PRESENT(u_minus_j)) THEN
727 5396 : IF (ASSOCIATED(qs_kind%dft_plus_u)) THEN
728 2698 : u_minus_j = qs_kind%dft_plus_u%u_minus_j
729 : ELSE
730 2698 : u_minus_j = 0.0_dp
731 : END IF
732 : END IF
733 77338383 : IF (PRESENT(u_minus_j_target)) THEN
734 5422 : IF (ASSOCIATED(qs_kind%dft_plus_u)) THEN
735 2698 : u_minus_j_target = qs_kind%dft_plus_u%u_minus_j_target
736 : ELSE
737 2724 : u_minus_j_target = 0.0_dp
738 : END IF
739 : END IF
740 77338383 : IF (PRESENT(U_of_dft_plus_u)) THEN
741 26 : IF (ASSOCIATED(qs_kind%dft_plus_u)) THEN
742 0 : U_of_dft_plus_u = qs_kind%dft_plus_u%U
743 : ELSE
744 26 : U_of_dft_plus_u = 0.0_dp
745 : END IF
746 : END IF
747 77338383 : IF (PRESENT(J_of_dft_plus_u)) THEN
748 26 : IF (ASSOCIATED(qs_kind%dft_plus_u)) THEN
749 0 : J_of_dft_plus_u = qs_kind%dft_plus_u%J
750 : ELSE
751 26 : J_of_dft_plus_u = 0.0_dp
752 : END IF
753 : END IF
754 77338383 : IF (PRESENT(alpha_of_dft_plus_u)) THEN
755 26 : IF (ASSOCIATED(qs_kind%dft_plus_u)) THEN
756 0 : alpha_of_dft_plus_u = qs_kind%dft_plus_u%alpha
757 : ELSE
758 26 : alpha_of_dft_plus_u = 0.0_dp
759 : END IF
760 : END IF
761 77338383 : IF (PRESENT(beta_of_dft_plus_u)) THEN
762 26 : IF (ASSOCIATED(qs_kind%dft_plus_u)) THEN
763 0 : beta_of_dft_plus_u = qs_kind%dft_plus_u%beta
764 : ELSE
765 26 : beta_of_dft_plus_u = 0.0_dp
766 : END IF
767 : END IF
768 77338383 : IF (PRESENT(J0_of_dft_plus_u)) THEN
769 26 : IF (ASSOCIATED(qs_kind%dft_plus_u)) THEN
770 0 : J0_of_dft_plus_u = qs_kind%dft_plus_u%J0
771 : ELSE
772 26 : J0_of_dft_plus_u = 0.0_dp
773 : END IF
774 : END IF
775 77338383 : IF (PRESENT(occupation_of_dft_plus_u)) THEN
776 26 : IF (ASSOCIATED(qs_kind%dft_plus_u)) THEN
777 0 : occupation_of_dft_plus_u = qs_kind%dft_plus_u%occupation
778 : ELSE
779 26 : occupation_of_dft_plus_u = -1.0_dp
780 : END IF
781 : END IF
782 :
783 77338383 : IF (PRESENT(init_u_ramping_each_scf)) THEN
784 184 : IF (ASSOCIATED(qs_kind%dft_plus_u)) THEN
785 92 : init_u_ramping_each_scf = qs_kind%dft_plus_u%init_u_ramping_each_scf
786 : ELSE
787 92 : init_u_ramping_each_scf = .FALSE.
788 : END IF
789 : END IF
790 77338383 : IF (PRESENT(u_ramping)) THEN
791 5580 : IF (ASSOCIATED(qs_kind%dft_plus_u)) THEN
792 2790 : u_ramping = qs_kind%dft_plus_u%u_ramping
793 : ELSE
794 2790 : u_ramping = 0.0_dp
795 : END IF
796 : END IF
797 77338383 : IF (PRESENT(eps_u_ramping)) THEN
798 5396 : IF (ASSOCIATED(qs_kind%dft_plus_u)) THEN
799 2698 : eps_u_ramping = qs_kind%dft_plus_u%eps_u_ramping
800 : ELSE
801 2698 : eps_u_ramping = 1.0E-5_dp
802 : END IF
803 : END IF
804 77338383 : IF (PRESENT(nelec)) THEN
805 3840 : NULLIFY (nelec)
806 3840 : IF (ASSOCIATED(qs_kind%dft_plus_u)) THEN
807 1920 : IF (ASSOCIATED(qs_kind%dft_plus_u%nelec)) THEN
808 0 : nelec => qs_kind%dft_plus_u%nelec
809 : END IF
810 : END IF
811 : END IF
812 77338383 : IF (PRESENT(orbitals)) THEN
813 4276 : NULLIFY (orbitals)
814 4276 : IF (ASSOCIATED(qs_kind%dft_plus_u)) THEN
815 2138 : IF (ASSOCIATED(qs_kind%dft_plus_u%orbitals)) THEN
816 136 : orbitals => qs_kind%dft_plus_u%orbitals
817 : END IF
818 : END IF
819 : END IF
820 77338383 : IF (PRESENT(eps_scf)) THEN
821 4276 : IF (ASSOCIATED(qs_kind%dft_plus_u)) THEN
822 2138 : eps_scf = qs_kind%dft_plus_u%eps_scf
823 : ELSE
824 2138 : eps_scf = 1.0E30_dp
825 : END IF
826 : END IF
827 77338383 : IF (PRESENT(max_scf)) THEN
828 4276 : IF (ASSOCIATED(qs_kind%dft_plus_u)) THEN
829 2138 : max_scf = qs_kind%dft_plus_u%max_scf
830 : ELSE
831 2138 : max_scf = -1
832 : END IF
833 : END IF
834 77338383 : IF (PRESENT(smear)) THEN
835 4276 : IF (ASSOCIATED(qs_kind%dft_plus_u)) THEN
836 2138 : smear = qs_kind%dft_plus_u%smear
837 : ELSE
838 2138 : smear = .FALSE.
839 : END IF
840 : END IF
841 77338383 : IF (PRESENT(dispersion)) dispersion => qs_kind%dispersion
842 77338383 : IF (PRESENT(bs_occupation)) bs_occupation = qs_kind%bs_occupation
843 77338383 : IF (PRESENT(addel)) addel => qs_kind%addel
844 77338383 : IF (PRESENT(laddel)) laddel => qs_kind%laddel
845 77338383 : IF (PRESENT(naddel)) naddel => qs_kind%naddel
846 :
847 77338383 : IF (PRESENT(magnetization)) magnetization = qs_kind%magnetization
848 :
849 77338383 : IF (PRESENT(no_optimize)) no_optimize = qs_kind%no_optimize
850 :
851 77338383 : IF (PRESENT(reltmat)) reltmat => qs_kind%reltmat
852 :
853 77338383 : IF (PRESENT(mao)) mao = qs_kind%mao
854 :
855 77338383 : IF (PRESENT(lmax_dftb)) lmax_dftb = qs_kind%lmax_dftb
856 :
857 77338383 : IF (PRESENT(pao_basis_size)) pao_basis_size = qs_kind%pao_basis_size
858 77338383 : IF (PRESENT(pao_model_file)) pao_model_file = qs_kind%pao_model_file
859 77338383 : IF (PRESENT(pao_potentials)) pao_potentials => qs_kind%pao_potentials
860 77338383 : IF (PRESENT(pao_descriptors)) pao_descriptors => qs_kind%pao_descriptors
861 77338383 : 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 4538219 : 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 4538219 : IF (PRESENT(basis_type)) THEN
948 4136637 : my_basis_type = basis_type
949 : ELSE
950 401582 : my_basis_type = "ORB"
951 : END IF
952 :
953 4538219 : IF (ASSOCIATED(qs_kind_set)) THEN
954 :
955 4538219 : IF (PRESENT(maxcgf)) maxcgf = 0
956 4538219 : IF (PRESENT(maxco)) maxco = 0
957 4538219 : IF (PRESENT(maxco_proj)) maxco_proj = 0
958 4538219 : IF (PRESENT(maxg_iso_not0)) maxg_iso_not0 = 0
959 4538219 : IF (PRESENT(maxgtops)) maxgtops = 0
960 4538219 : IF (PRESENT(maxlgto)) maxlgto = -1
961 4538219 : IF (PRESENT(maxlppl)) maxlppl = -1
962 4538219 : IF (PRESENT(maxlppnl)) maxlppnl = -1
963 4538219 : IF (PRESENT(maxpol)) maxpol = -1
964 4538219 : IF (PRESENT(maxlprj)) maxlprj = -1
965 4538219 : IF (PRESENT(maxnset)) maxnset = 0
966 4538219 : IF (PRESENT(maxppnl)) maxppnl = 0
967 4538219 : IF (PRESENT(maxsgf)) maxsgf = 0
968 4538219 : IF (PRESENT(maxsgf_set)) maxsgf_set = 0
969 4538219 : IF (PRESENT(ncgf)) ncgf = 0
970 4538219 : IF (PRESENT(nelectron)) nelectron = 0
971 4538219 : IF (PRESENT(npgf)) npgf = 0
972 4538219 : IF (PRESENT(nset)) nset = 0
973 4538219 : IF (PRESENT(nsgf)) nsgf = 0
974 4538219 : IF (PRESENT(nshell)) nshell = 0
975 4538219 : IF (PRESENT(all_potential_present)) all_potential_present = .FALSE.
976 4538219 : IF (PRESENT(tnadd_potential_present)) tnadd_potential_present = .FALSE.
977 4538219 : IF (PRESENT(gth_potential_present)) gth_potential_present = .FALSE.
978 4538219 : IF (PRESENT(sgp_potential_present)) sgp_potential_present = .FALSE.
979 4538219 : IF (PRESENT(cneo_potential_present)) cneo_potential_present = .FALSE.
980 4538219 : IF (PRESENT(nkind_q)) nkind_q = 0
981 4538219 : IF (PRESENT(natom_q)) natom_q = 0
982 4538219 : IF (PRESENT(paw_atom_present)) paw_atom_present = .FALSE.
983 4538219 : IF (PRESENT(max_ngrid_rad)) max_ngrid_rad = 0
984 4538219 : IF (PRESENT(max_sph_harm)) max_sph_harm = 0
985 4538219 : IF (PRESENT(lmax_rho0)) lmax_rho0 = 0
986 4538219 : IF (PRESENT(basis_rcut)) basis_rcut = 0.0_dp
987 4538219 : IF (PRESENT(total_zeff_corr)) total_zeff_corr = 0.0_dp
988 4538219 : IF (PRESENT(npgf_seg)) npgf_seg = 0
989 :
990 4538219 : nkind = SIZE(qs_kind_set)
991 13711333 : DO ikind = 1, nkind
992 9173114 : 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 9173114 : lmax_rho0=lmax_rho0_kind)
1007 :
1008 9173114 : IF (PRESENT(maxlppl) .AND. ASSOCIATED(gth_potential)) THEN
1009 48888 : CALL get_potential(potential=gth_potential, nexp_ppl=n)
1010 48888 : maxlppl = MAX(maxlppl, 2*(n - 1))
1011 10891 : ELSE IF (PRESENT(maxlppl) .AND. ASSOCIATED(sgp_potential)) THEN
1012 324 : CALL get_potential(potential=sgp_potential, nrloc=nrloc, ecp_semi_local=ecp_semi_local)
1013 3564 : n = MAXVAL(nrloc) - 2
1014 324 : maxlppl = MAX(maxlppl, 2*(n - 1))
1015 324 : IF (ecp_semi_local) THEN
1016 294 : CALL get_potential(potential=sgp_potential, sl_lmax=imax, nrpot=nrpot)
1017 52038 : n = MAXVAL(nrpot) - 2
1018 294 : n = 2*(n - 1) + imax
1019 294 : maxlppl = MAX(maxlppl, n)
1020 : END IF
1021 : END IF
1022 :
1023 9173114 : IF (PRESENT(maxlppnl) .AND. ASSOCIATED(gth_potential)) THEN
1024 45675 : CALL get_potential(potential=gth_potential, lprj_ppnl_max=imax)
1025 45675 : maxlppnl = MAX(maxlppnl, imax)
1026 10687 : ELSE IF (PRESENT(maxlppnl) .AND. ASSOCIATED(sgp_potential)) THEN
1027 142 : CALL get_potential(potential=sgp_potential, lmax=imax)
1028 142 : maxlppnl = MAX(maxlppnl, imax)
1029 : END IF
1030 :
1031 9173114 : IF (PRESENT(maxpol) .AND. ASSOCIATED(tnadd_potential)) THEN
1032 78 : CALL get_potential(potential=tnadd_potential, npol=n)
1033 78 : maxpol = MAX(maxpol, 2*(n - 1))
1034 : END IF
1035 :
1036 9173114 : IF (PRESENT(maxco_proj) .AND. ASSOCIATED(paw_proj_set)) THEN
1037 6594 : CALL get_paw_proj_set(paw_proj_set=paw_proj_set, ncgauprj=imax)
1038 6594 : maxco_proj = MAX(maxco_proj, imax)
1039 : END IF
1040 :
1041 9173114 : IF (PRESENT(maxlprj) .AND. ASSOCIATED(paw_proj_set)) THEN
1042 6594 : CALL get_paw_proj_set(paw_proj_set=paw_proj_set, maxl=imax)
1043 6594 : maxlprj = MAX(maxlprj, imax)
1044 : END IF
1045 :
1046 9173114 : IF (PRESENT(maxppnl) .AND. ASSOCIATED(gth_potential)) THEN
1047 31422 : CALL get_potential(potential=gth_potential, nppnl=imax)
1048 31422 : maxppnl = MAX(maxppnl, imax)
1049 248 : ELSE IF (PRESENT(maxppnl) .AND. ASSOCIATED(sgp_potential)) THEN
1050 10 : CALL get_potential(potential=sgp_potential, nppnl=imax)
1051 10 : maxppnl = MAX(maxppnl, imax)
1052 : END IF
1053 :
1054 9173114 : IF (my_basis_type(1:3) == "NUC" .AND. .NOT. ASSOCIATED(cneo_potential)) THEN
1055 8305 : NULLIFY (tmp_basis_set)
1056 : ELSE
1057 : CALL get_basis_from_container(qs_kind%basis_sets, basis_set=tmp_basis_set, &
1058 9164809 : basis_type=my_basis_type)
1059 : END IF
1060 :
1061 9173114 : IF (PRESENT(maxcgf)) THEN
1062 0 : IF (ASSOCIATED(tmp_basis_set)) THEN
1063 0 : CALL get_gto_basis_set(gto_basis_set=tmp_basis_set, ncgf=imax)
1064 0 : maxcgf = MAX(maxcgf, imax)
1065 0 : ELSE IF (ASSOCIATED(qs_kind%dftb_parameter)) THEN
1066 0 : CALL get_dftb_atom_param(dftb_parameter=dftb_parameter, lmax=imax)
1067 0 : imax = ((imax + 1)*(imax + 2)*(imax + 3))/6
1068 0 : maxcgf = MAX(maxcgf, imax)
1069 : END IF
1070 : END IF
1071 :
1072 9173114 : IF (PRESENT(maxco)) THEN
1073 8200050 : IF (ASSOCIATED(tmp_basis_set)) THEN
1074 8199942 : IF (PRESENT(maxder)) THEN
1075 : CALL get_gto_basis_set(gto_basis_set=tmp_basis_set, &
1076 0 : maxco=imax, maxder=maxder)
1077 : ELSE
1078 8199942 : CALL get_gto_basis_set(gto_basis_set=tmp_basis_set, maxco=imax)
1079 : END IF
1080 8199942 : maxco = MAX(maxco, imax)
1081 : END IF
1082 8200050 : IF (ASSOCIATED(gth_potential)) THEN
1083 726002 : CALL get_potential(potential=gth_potential, lprj_ppnl_max=imax)
1084 726002 : maxco = MAX(maxco, ncoset(imax))
1085 : END IF
1086 8200050 : IF (ASSOCIATED(sgp_potential)) THEN
1087 2366 : CALL get_potential(potential=sgp_potential, lmax=imax)
1088 2366 : maxco = MAX(maxco, ncoset(imax))
1089 2366 : CALL get_potential(potential=sgp_potential, sl_lmax=imax)
1090 2366 : maxco = MAX(maxco, ncoset(imax))
1091 : END IF
1092 : END IF
1093 :
1094 9173114 : IF (PRESENT(maxgtops)) THEN
1095 161126 : IF (ASSOCIATED(tmp_basis_set)) THEN
1096 161126 : CALL get_gto_basis_set(gto_basis_set=tmp_basis_set, maxso=imax, nset=n)
1097 161126 : maxgtops = MAX(maxgtops, n*imax)
1098 : END IF
1099 : END IF
1100 :
1101 9173114 : IF (PRESENT(maxlgto)) THEN
1102 7758517 : IF (ASSOCIATED(tmp_basis_set)) THEN
1103 7721965 : CALL get_gto_basis_set(gto_basis_set=tmp_basis_set, maxl=imax)
1104 7721965 : maxlgto = MAX(maxlgto, imax)
1105 36552 : ELSE IF (ASSOCIATED(qs_kind%dftb_parameter)) THEN
1106 5483 : CALL get_dftb_atom_param(dftb_parameter=dftb_parameter, lmax=imax)
1107 5483 : maxlgto = MAX(maxlgto, imax)
1108 : END IF
1109 : END IF
1110 :
1111 9173114 : IF (PRESENT(maxnset)) THEN
1112 80359 : IF (ASSOCIATED(tmp_basis_set)) THEN
1113 80347 : CALL get_gto_basis_set(gto_basis_set=tmp_basis_set, nset=n)
1114 80347 : maxnset = MAX(maxnset, n)
1115 : END IF
1116 : END IF
1117 :
1118 9173114 : IF (PRESENT(maxsgf)) THEN
1119 8019042 : IF (ASSOCIATED(tmp_basis_set)) THEN
1120 8019006 : CALL get_gto_basis_set(gto_basis_set=tmp_basis_set, nsgf=imax)
1121 8019006 : maxsgf = MAX(maxsgf, imax)
1122 : END IF
1123 : END IF
1124 :
1125 9173114 : IF (PRESENT(maxsgf_set)) THEN
1126 502062 : IF (ASSOCIATED(tmp_basis_set)) THEN
1127 501974 : CALL get_gto_basis_set(gto_basis_set=tmp_basis_set, maxsgf_set=imax)
1128 501974 : maxsgf_set = MAX(maxsgf_set, imax)
1129 : END IF
1130 : END IF
1131 :
1132 9173114 : IF (PRESENT(ncgf)) THEN
1133 45926 : IF (ASSOCIATED(tmp_basis_set)) THEN
1134 13317 : CALL get_gto_basis_set(gto_basis_set=tmp_basis_set, ncgf=n)
1135 13317 : ncgf = ncgf + n*qs_kind_set(ikind)%natom
1136 32609 : ELSE IF (ASSOCIATED(qs_kind%dftb_parameter)) THEN
1137 1572 : CALL get_dftb_atom_param(dftb_parameter=dftb_parameter, lmax=imax)
1138 1572 : n = ((imax + 1)*(imax + 2)*(imax + 3))/6
1139 1572 : ncgf = ncgf + n*qs_kind_set(ikind)%natom
1140 : END IF
1141 : END IF
1142 :
1143 9173114 : IF (PRESENT(npgf)) THEN
1144 41015 : IF (ASSOCIATED(tmp_basis_set)) THEN
1145 8433 : CALL get_gto_basis_set(gto_basis_set=tmp_basis_set, npgf_sum=n)
1146 8433 : npgf = npgf + n*qs_kind_set(ikind)%natom
1147 : END IF
1148 : END IF
1149 :
1150 9173114 : IF (PRESENT(nset)) THEN
1151 41015 : IF (ASSOCIATED(tmp_basis_set)) THEN
1152 8433 : CALL get_gto_basis_set(gto_basis_set=tmp_basis_set, nset=n)
1153 8433 : nset = nset + n*qs_kind_set(ikind)%natom
1154 : END IF
1155 : END IF
1156 :
1157 9173114 : IF (PRESENT(nsgf)) THEN
1158 124751 : IF (ASSOCIATED(tmp_basis_set)) THEN
1159 72526 : CALL get_gto_basis_set(gto_basis_set=tmp_basis_set, nsgf=n)
1160 72526 : nsgf = nsgf + n*qs_kind_set(ikind)%natom
1161 52225 : ELSE IF (ASSOCIATED(qs_kind%dftb_parameter)) THEN
1162 21186 : CALL get_dftb_atom_param(dftb_parameter=dftb_parameter, natorb=n)
1163 21186 : nsgf = nsgf + n*qs_kind_set(ikind)%natom
1164 : END IF
1165 : END IF
1166 :
1167 9173114 : IF (PRESENT(nshell)) THEN
1168 41037 : IF (ASSOCIATED(tmp_basis_set)) THEN
1169 8455 : CALL get_gto_basis_set(gto_basis_set=tmp_basis_set, nshell_sum=n)
1170 8455 : nshell = nshell + n*qs_kind_set(ikind)%natom
1171 32582 : ELSE IF (ASSOCIATED(qs_kind%dftb_parameter)) THEN
1172 1545 : CALL get_dftb_atom_param(dftb_parameter=dftb_parameter, lmax=n)
1173 1545 : nshell = nshell + (n + 1)*qs_kind_set(ikind)%natom
1174 : END IF
1175 : END IF
1176 :
1177 9173114 : IF (PRESENT(nelectron)) THEN
1178 242334 : IF (ASSOCIATED(qs_kind%all_potential)) THEN
1179 : CALL get_potential(potential=qs_kind%all_potential, &
1180 23726 : zeff=zeff, zeff_correction=zeff_correction)
1181 218608 : ELSE IF (ASSOCIATED(qs_kind%gth_potential)) THEN
1182 : CALL get_potential(potential=qs_kind%gth_potential, &
1183 216558 : zeff=zeff, zeff_correction=zeff_correction)
1184 2050 : ELSE IF (ASSOCIATED(qs_kind%sgp_potential)) THEN
1185 : CALL get_potential(potential=qs_kind%sgp_potential, &
1186 1100 : zeff=zeff, zeff_correction=zeff_correction)
1187 950 : ELSE IF (ASSOCIATED(qs_kind%cneo_potential)) THEN
1188 : CALL get_cneo_potential(potential=qs_kind%cneo_potential, &
1189 56 : zeff=zeff)
1190 56 : zeff_correction = 0.0_dp
1191 : ELSE
1192 894 : zeff = 0.0_dp
1193 894 : zeff_correction = 0.0_dp
1194 : END IF
1195 242334 : nelectron = nelectron + qs_kind_set(ikind)%natom*NINT(zeff - zeff_correction)
1196 : END IF
1197 :
1198 9173114 : IF (PRESENT(basis_rcut)) THEN
1199 286 : IF (ASSOCIATED(tmp_basis_set)) THEN
1200 0 : CALL get_gto_basis_set(gto_basis_set=tmp_basis_set, kind_radius=brcut)
1201 0 : basis_rcut = MAX(basis_rcut, brcut)
1202 286 : ELSE IF (ASSOCIATED(qs_kind%dftb_parameter)) THEN
1203 286 : CALL get_dftb_atom_param(dftb_parameter=dftb_parameter, cutoff=brcut)
1204 286 : basis_rcut = MAX(basis_rcut, brcut)
1205 : END IF
1206 : END IF
1207 :
1208 9173114 : IF (PRESENT(total_zeff_corr)) THEN
1209 16365 : IF (ASSOCIATED(qs_kind%all_potential)) THEN
1210 : CALL get_potential(potential=qs_kind%all_potential, &
1211 6712 : zeff=zeff, zeff_correction=zeff_correction)
1212 9653 : ELSE IF (ASSOCIATED(qs_kind%gth_potential)) THEN
1213 : CALL get_potential(potential=qs_kind%gth_potential, &
1214 9405 : zeff=zeff, zeff_correction=zeff_correction)
1215 248 : ELSE IF (ASSOCIATED(qs_kind%sgp_potential)) THEN
1216 : CALL get_potential(potential=qs_kind%sgp_potential, &
1217 88 : zeff=zeff, zeff_correction=zeff_correction)
1218 : ELSE
1219 160 : zeff = 0.0_dp
1220 160 : zeff_correction = 0.0_dp
1221 : END IF
1222 16365 : total_zeff_corr = total_zeff_corr + qs_kind_set(ikind)%natom*zeff_correction
1223 : END IF
1224 :
1225 9173114 : IF (PRESENT(all_potential_present)) THEN
1226 79885 : IF (ASSOCIATED(all_potential)) THEN
1227 49088 : all_potential_present = .TRUE.
1228 : END IF
1229 : END IF
1230 :
1231 9173114 : IF (PRESENT(tnadd_potential_present)) THEN
1232 0 : IF (ASSOCIATED(tnadd_potential)) THEN
1233 0 : tnadd_potential_present = .TRUE.
1234 : END IF
1235 : END IF
1236 :
1237 9173114 : IF (PRESENT(gth_potential_present)) THEN
1238 63464 : IF (ASSOCIATED(gth_potential)) THEN
1239 20812 : gth_potential_present = .TRUE.
1240 : END IF
1241 : END IF
1242 :
1243 9173114 : IF (PRESENT(sgp_potential_present)) THEN
1244 63470 : IF (ASSOCIATED(sgp_potential)) THEN
1245 138 : sgp_potential_present = .TRUE.
1246 : END IF
1247 : END IF
1248 :
1249 9173114 : IF (PRESENT(cneo_potential_present)) THEN
1250 86599 : IF (ASSOCIATED(cneo_potential)) THEN
1251 24 : cneo_potential_present = .TRUE.
1252 : END IF
1253 : END IF
1254 :
1255 9173114 : IF (PRESENT(nkind_q)) THEN
1256 8203 : IF (ASSOCIATED(cneo_potential)) THEN
1257 4 : nkind_q = nkind_q + 1
1258 : END IF
1259 : END IF
1260 :
1261 9173114 : IF (PRESENT(natom_q)) THEN
1262 8203 : IF (ASSOCIATED(cneo_potential)) THEN
1263 4 : natom_q = natom_q + qs_kind_set(ikind)%natom
1264 : END IF
1265 : END IF
1266 :
1267 9173114 : IF (PRESENT(paw_atom_present)) THEN
1268 63464 : IF (paw_atom) THEN
1269 4656 : paw_atom_present = .TRUE.
1270 : END IF
1271 : END IF
1272 :
1273 9173114 : IF (PRESENT(dft_plus_u_atom_present)) THEN
1274 16365 : IF (dft_plus_u_atom) THEN
1275 36 : dft_plus_u_atom_present = .TRUE.
1276 : END IF
1277 : END IF
1278 :
1279 9173114 : IF (PRESENT(max_ngrid_rad)) THEN
1280 0 : max_ngrid_rad = MAX(max_ngrid_rad, ngrid_rad)
1281 : END IF
1282 :
1283 9173114 : IF (PRESENT(max_sph_harm)) THEN
1284 0 : max_sph_harm = MAX(max_sph_harm, max_s_harm)
1285 : END IF
1286 :
1287 9173114 : IF (PRESENT(maxg_iso_not0)) THEN
1288 50500 : maxg_iso_not0 = MAX(maxg_iso_not0, max_iso_not0)
1289 : END IF
1290 :
1291 9173114 : IF (PRESENT(lmax_rho0)) THEN
1292 0 : lmax_rho0 = MAX(lmax_rho0, lmax_rho0_kind)
1293 : END IF
1294 :
1295 22884447 : IF (PRESENT(npgf_seg)) THEN
1296 22 : IF (ASSOCIATED(tmp_basis_set)) THEN
1297 22 : CALL get_gto_basis_set(gto_basis_set=tmp_basis_set, npgf_seg_sum=n)
1298 22 : npgf_seg = npgf_seg + n*qs_kind_set(ikind)%natom
1299 : END IF
1300 : END IF
1301 :
1302 : END DO
1303 : ELSE
1304 0 : CPABORT("The pointer qs_kind_set is not associated")
1305 : END IF
1306 :
1307 4538219 : 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 16365 : 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 16365 : CALL timeset(routineN, handle)
1325 :
1326 16365 : CPASSERT(ASSOCIATED(qs_kind))
1327 :
1328 16365 : IF (ASSOCIATED(qs_kind%gth_potential)) THEN
1329 9405 : CALL init_potential(qs_kind%gth_potential)
1330 6960 : ELSE IF (ASSOCIATED(qs_kind%sgp_potential)) THEN
1331 88 : CALL init_potential(qs_kind%sgp_potential)
1332 : END IF
1333 :
1334 343665 : DO i = 1, SIZE(qs_kind%basis_sets, 1)
1335 327300 : NULLIFY (tmp_basis_set)
1336 : CALL get_basis_from_container(qs_kind%basis_sets, basis_set=tmp_basis_set, &
1337 327300 : inumbas=i, basis_type=basis_type)
1338 327300 : IF (basis_type == "") CYCLE
1339 34755 : IF (basis_type == "AUX") THEN
1340 0 : IF (tmp_basis_set%norm_type < 0) tmp_basis_set%norm_type = 1
1341 0 : CALL init_aux_basis_set(tmp_basis_set)
1342 : ELSE
1343 18390 : IF (tmp_basis_set%norm_type < 0) tmp_basis_set%norm_type = 2
1344 18390 : CALL init_orb_basis_set(tmp_basis_set)
1345 : END IF
1346 : END DO
1347 :
1348 16365 : CALL timestop(handle)
1349 :
1350 16365 : 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 8664 : 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 8664 : CALL timeset(routineN, handle)
1368 :
1369 8664 : IF (.NOT. ASSOCIATED(qs_kind_set)) THEN
1370 0 : CPABORT("init_qs_kind_set: The pointer qs_kind_set is not associated")
1371 : END IF
1372 :
1373 25029 : DO ikind = 1, SIZE(qs_kind_set)
1374 16365 : qs_kind => qs_kind_set(ikind)
1375 25029 : CALL init_qs_kind(qs_kind)
1376 : END DO
1377 :
1378 8664 : CALL timestop(handle)
1379 :
1380 8664 : 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 1502 : 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 1502 : my_mod_control = .TRUE.
1405 1502 : IF (PRESENT(modify_qs_control)) THEN
1406 146 : my_mod_control = modify_qs_control
1407 : END IF
1408 :
1409 1502 : IF (ASSOCIATED(qs_kind_set)) THEN
1410 :
1411 1502 : IF (my_mod_control) qs_control%gapw_control%non_paw_atoms = .FALSE.
1412 1502 : nkind = SIZE(qs_kind_set)
1413 :
1414 4290 : DO ikind = 1, nkind
1415 :
1416 2788 : qs_kind => qs_kind_set(ikind)
1417 :
1418 2788 : 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 2788 : max_rad_local=max_rad_local_type, gpw_type_forced=gpw)
1421 :
1422 2788 : NULLIFY (soft_basis)
1423 2788 : 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 5444 : ASSOCIATED(qs_kind%cneo_potential)), gpw)
1431 2788 : CALL add_basis_set_to_container(qs_kind%basis_sets, soft_basis, "ORB_SOFT")
1432 2788 : CALL set_qs_kind(qs_kind=qs_kind, paw_atom=paw_atom)
1433 :
1434 2788 : bas1c = qs_control%gapw_control%basis_1c
1435 2788 : NULLIFY (basis_1c)
1436 2564 : SELECT CASE (bas1c)
1437 : CASE (gapw_1c_orb)
1438 2564 : ilevel = 0
1439 : CASE (gapw_1c_small)
1440 208 : ilevel = 1
1441 : CASE (gapw_1c_medium)
1442 4 : ilevel = 2
1443 : CASE (gapw_1c_large)
1444 8 : ilevel = 3
1445 : CASE (gapw_1c_very_large)
1446 4 : ilevel = 4
1447 : CASE DEFAULT
1448 2788 : CPABORT("basis_1c type")
1449 : END SELECT
1450 2788 : CALL remove_basis_from_container(qs_kind%basis_sets, basis_type="GAPW_1C")
1451 2788 : CALL create_1c_basis(orb_basis, soft_basis, basis_1c, ilevel)
1452 2788 : CALL get_gto_basis_set(gto_basis_set=orb_basis, name=bsname)
1453 2788 : basis_1c%name = TRIM(bsname)//"_1c"
1454 2788 : CALL add_basis_set_to_container(qs_kind%basis_sets, basis_1c, "GAPW_1C")
1455 2788 : IF (paw_atom) THEN
1456 2376 : CALL allocate_paw_proj_set(qs_kind%paw_proj_set)
1457 2376 : 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 2376 : max_rad_local_type, force_env_section)
1460 : ELSE
1461 412 : 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 2788 : NULLIFY (qs_kind%grid_atom, qs_kind%harmonics)
1466 2788 : CALL allocate_grid_atom(qs_kind%grid_atom)
1467 7078 : CALL allocate_harmonics_atom(qs_kind%harmonics)
1468 :
1469 : END DO
1470 :
1471 1502 : IF (my_mod_control) THEN
1472 1356 : IF (qs_control%gapw_control%non_paw_atoms) THEN
1473 202 : qs_control%gapw_control%nopaw_as_gpw = .TRUE.
1474 : ELSE
1475 1154 : qs_control%gapw_control%nopaw_as_gpw = .FALSE.
1476 : END IF
1477 : END IF
1478 : ELSE
1479 0 : CPABORT("The pointer qs_kind_set is not associated")
1480 : END IF
1481 :
1482 1502 : END SUBROUTINE init_gapw_basis_set
1483 : ! **************************************************************************************************
1484 : !> \brief ...
1485 : !> \param qs_kind_set ...
1486 : ! **************************************************************************************************
1487 1502 : 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 1502 : INTEGER, DIMENSION(:), POINTER :: nct_nlcc
1494 : LOGICAL :: nlcc, nlcc_type, paw_atom
1495 : REAL(dp) :: alpha, coa, cval
1496 1502 : REAL(KIND=dp), DIMENSION(:), POINTER :: a_nlcc, alpha_nlcc, c_nlcc, fe, rc, rr
1497 1502 : 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 1502 : IF (ASSOCIATED(qs_kind_set)) THEN
1503 1502 : nlcc = has_nlcc(qs_kind_set)
1504 1502 : IF (nlcc) THEN
1505 10 : nkind = SIZE(qs_kind_set)
1506 28 : DO ikind = 1, nkind
1507 18 : qs_kind => qs_kind_set(ikind)
1508 18 : CALL get_qs_kind(qs_kind, paw_atom=paw_atom)
1509 28 : IF (paw_atom) THEN
1510 16 : CALL get_qs_kind(qs_kind, gth_potential=gth_potential)
1511 16 : CALL get_qs_kind(qs_kind, sgp_potential=sgp_potential)
1512 16 : IF (ASSOCIATED(gth_potential)) THEN
1513 : CALL get_potential(potential=gth_potential, nlcc_present=nlcc_type, &
1514 16 : nexp_nlcc=nexp_nlcc, alpha_nlcc=alpha_nlcc, nct_nlcc=nct_nlcc, cval_nlcc=cval_nlcc)
1515 16 : IF (nlcc_type) THEN
1516 10 : nr = qs_kind%grid_atom%nr
1517 10 : rr => qs_kind%grid_atom%rad
1518 60 : ALLOCATE (qs_kind%nlcc_pot(nr, 2), rc(nr), fe(nr))
1519 30 : den => qs_kind%nlcc_pot
1520 1030 : den = 0.0_dp
1521 20 : DO i = 1, nexp_nlcc
1522 10 : alpha = alpha_nlcc(i)
1523 1010 : rc(:) = rr(:)/alpha
1524 1010 : fe(:) = EXP(-0.5_dp*rc(:)*rc(:))
1525 10 : nc = nct_nlcc(i)
1526 32 : DO ic = 1, nc
1527 12 : cval = cval_nlcc(ic, i)
1528 12 : coa = cval/alpha
1529 1212 : den(:, 1) = den(:, 1) + fe(:)*rc**(2*ic - 2)*cval
1530 1212 : den(:, 2) = den(:, 2) - fe(:)*rc**(2*ic - 1)*coa
1531 22 : IF (ic > 1) THEN
1532 202 : 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 10 : DEALLOCATE (rc, fe)
1537 : END IF
1538 0 : ELSE IF (ASSOCIATED(sgp_potential)) THEN
1539 : CALL get_potential(potential=sgp_potential, has_nlcc=nlcc_type, &
1540 0 : n_nlcc=n_nlcc, a_nlcc=a_nlcc, c_nlcc=c_nlcc)
1541 0 : IF (nlcc_type) THEN
1542 0 : nr = qs_kind%grid_atom%nr
1543 0 : rr => qs_kind%grid_atom%rad
1544 0 : ALLOCATE (qs_kind%nlcc_pot(nr, 2), rc(nr), fe(nr))
1545 0 : den => qs_kind%nlcc_pot
1546 0 : den = 0.0_dp
1547 0 : DO i = 1, n_nlcc
1548 0 : alpha = a_nlcc(i)
1549 0 : fe(:) = EXP(-alpha*rr(:)*rr(:))
1550 0 : cval = c_nlcc(i)
1551 0 : den(:, 1) = den(:, 1) + cval*fe(:)
1552 0 : den(:, 2) = den(:, 2) - 2.0_dp*alpha*cval*rr(:)*fe(:)
1553 : END DO
1554 0 : 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 0 : CPABORT("The pointer qs_kind_set is not associated")
1564 : END IF
1565 :
1566 1502 : END SUBROUTINE init_gapw_nlcc
1567 :
1568 : ! **************************************************************************************************
1569 : !> \brief ...
1570 : !> \param qs_kind_set ...
1571 : !> \param qs_control ...
1572 : ! **************************************************************************************************
1573 8 : 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 8 : IF (ASSOCIATED(qs_kind_set)) THEN
1585 :
1586 8 : nkind = SIZE(qs_kind_set)
1587 :
1588 22 : DO ikind = 1, nkind
1589 :
1590 14 : qs_kind => qs_kind_set(ikind)
1591 :
1592 22 : IF (ASSOCIATED(qs_kind%cneo_potential)) THEN
1593 : CALL get_qs_kind(qs_kind=qs_kind, basis_set=orb_basis, basis_type="NUC", &
1594 8 : hard_radius=rc)
1595 :
1596 8 : NULLIFY (soft_basis)
1597 8 : 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 8 : rc, paw_atom, .TRUE., .FALSE.)
1602 8 : 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 0 : CPABORT("The pointer qs_kind_set is not associated")
1609 : END IF
1610 :
1611 8 : 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 16447 : 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 16447 : 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 32894 : 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 32894 : REAL(KIND=dp), DIMENSION(:), POINTER :: a_nl, aloc, anlcc, cloc, cnlcc, nelec
1663 16447 : REAL(KIND=dp), DIMENSION(:, :), POINTER :: h_nl
1664 16447 : REAL(KIND=dp), DIMENSION(:, :, :), POINTER :: c_nl
1665 : TYPE(atom_ecppot_type) :: ecppot
1666 : TYPE(atom_sgp_potential_type) :: sgppot
1667 1726935 : 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 16447 : CALL timeset(routineN, handle)
1677 :
1678 16447 : NULLIFY (logger)
1679 16447 : logger => cp_get_default_logger()
1680 16447 : iounit = cp_logger_get_default_io_unit(logger)
1681 :
1682 16447 : NULLIFY (elec_conf)
1683 :
1684 16447 : update_input = .TRUE.
1685 345387 : basis_set_name(:) = ""
1686 345387 : basis_set_type(:) = ""
1687 345387 : basis_set_form(:) = ""
1688 16447 : potential_name = ""
1689 16447 : potential_type = ""
1690 16447 : kgpot_name = ""
1691 16447 : kgpot_type = ""
1692 16447 : z = -1
1693 16447 : zeff_correction = 0.0_dp
1694 16447 : explicit = .FALSE.
1695 16447 : explicit_basis = .FALSE.
1696 16447 : explicit_J = .FALSE.
1697 16447 : explicit_kgpot = .FALSE.
1698 16447 : explicit_potential = .FALSE.
1699 16447 : explicit_U = .FALSE.
1700 16447 : explicit_u_m_j = .FALSE.
1701 :
1702 16447 : dft_section => section_vals_get_subs_vals(force_env_section, "DFT")
1703 16447 : CALL section_vals_get(kind_section, n_repetition=n_rep)
1704 16447 : k_rep = -1
1705 16447 : akind_name = qs_kind%name
1706 16447 : CALL uppercase(akind_name)
1707 : ! First we use the atom_name to find out the proper KIND section
1708 23638 : DO i_rep = 1, n_rep
1709 : CALL section_vals_val_get(kind_section, "_SECTION_PARAMETERS_", &
1710 18296 : c_val=keyword, i_rep_section=i_rep)
1711 18296 : CALL uppercase(keyword)
1712 23638 : IF (keyword == akind_name) THEN
1713 11105 : k_rep = i_rep
1714 11105 : EXIT
1715 : END IF
1716 : END DO
1717 : ! The search for the KIND section failed.. check for a QM/MM link atom
1718 16447 : IF (k_rep < 1) THEN
1719 5342 : ipos = INDEX(qs_kind%name, "_")
1720 5342 : 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 64 : akind_name = qs_kind%name(1:ipos - 1)
1724 64 : CALL uppercase(akind_name)
1725 64 : DO i_rep = 1, n_rep
1726 : CALL section_vals_val_get(kind_section, "_SECTION_PARAMETERS_", &
1727 52 : c_val=keyword, i_rep_section=i_rep)
1728 52 : CALL uppercase(keyword)
1729 64 : IF (keyword == akind_name) THEN
1730 52 : k_rep = i_rep
1731 52 : 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 16447 : 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 5290 : element_symbol = qs_kind%element_symbol(1:2)
1741 5290 : CALL uppercase(element_symbol)
1742 5384 : DO i_rep = 1, n_rep
1743 : CALL section_vals_val_get(kind_section, "_SECTION_PARAMETERS_", &
1744 118 : c_val=keyword, i_rep_section=i_rep)
1745 118 : CALL uppercase(keyword)
1746 5384 : IF (keyword == element_symbol) THEN
1747 24 : k_rep = i_rep
1748 24 : 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 16447 : IF (k_rep < 1) THEN
1755 5282 : DO i_rep = 1, n_rep
1756 : CALL section_vals_val_get(kind_section, "_SECTION_PARAMETERS_", &
1757 72 : c_val=keyword, i_rep_section=i_rep)
1758 72 : CALL uppercase(keyword)
1759 5282 : IF (keyword == "DEFAULT") THEN
1760 56 : update_input = .FALSE.
1761 56 : k_rep = i_rep
1762 56 : EXIT
1763 : END IF
1764 : END DO
1765 : END IF
1766 16447 : 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 0 : " and for the DEFAULT section but no match was found. Check your input file!")
1772 : END IF
1773 : ! Retrieve information on element
1774 16447 : CALL get_ptable_info(qs_kind%element_symbol, ielement=z)
1775 :
1776 : ! Normal parsing of the KIND section
1777 16447 : 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 11237 : n_rep_val=nb_rep)
1783 11237 : IF (.NOT. explicit) nb_rep = 0
1784 11237 : CPASSERT(nb_rep <= maxbas)
1785 24447 : DO i = 1, nb_rep
1786 : CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
1787 13210 : keyword_name="BASIS_SET", i_rep_val=i, c_vals=tmpstringlist)
1788 13210 : IF (SIZE(tmpstringlist) == 1) THEN
1789 : ! default is orbital type and GTO
1790 9829 : basis_set_type(i) = "ORB"
1791 9829 : basis_set_form(i) = "GTO"
1792 9829 : basis_set_name(i) = tmpstringlist(1)
1793 3381 : ELSE IF (SIZE(tmpstringlist) == 2) THEN
1794 : ! default is GTO
1795 3377 : basis_set_type(i) = tmpstringlist(1)
1796 3377 : basis_set_form(i) = "GTO"
1797 3377 : basis_set_name(i) = tmpstringlist(2)
1798 4 : ELSE IF (SIZE(tmpstringlist) == 3) THEN
1799 4 : basis_set_type(i) = tmpstringlist(1)
1800 4 : basis_set_form(i) = tmpstringlist(2)
1801 4 : basis_set_name(i) = tmpstringlist(3)
1802 : ELSE
1803 : CALL cp_abort(__LOCATION__, &
1804 0 : "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 24447 : IF (basis_set_form(i) /= "GTO" .AND. basis_set_form(i) /= "STO") THEN
1808 0 : 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 11237 : 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 11237 : explicit=explicit)
1817 11237 : IF (explicit) THEN
1818 : CALL section_vals_val_get(kind_section, keyword_name="PAO_MODEL_FILE", i_rep_section=k_rep, &
1819 8 : c_val=qs_kind%pao_model_file)
1820 : END IF
1821 :
1822 : ! parse PAO_POTENTIAL sections
1823 11237 : pao_pot_section => section_vals_get_subs_vals(kind_section, "PAO_POTENTIAL", i_rep_section=k_rep)
1824 11237 : CALL section_vals_get(pao_pot_section, n_repetition=npaopot)
1825 22596 : ALLOCATE (qs_kind%pao_potentials(npaopot))
1826 11299 : DO ipaopot = 1, npaopot
1827 : CALL section_vals_val_get(pao_pot_section, keyword_name="MAXL", i_rep_section=ipaopot, &
1828 62 : 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 62 : 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 62 : 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 11299 : r_val=qs_kind%pao_potentials(ipaopot)%weight)
1835 : END DO
1836 :
1837 : ! parse PAO_DESCRIPTOR sections
1838 11237 : pao_desc_section => section_vals_get_subs_vals(kind_section, "PAO_DESCRIPTOR", i_rep_section=k_rep)
1839 11237 : CALL section_vals_get(pao_desc_section, n_repetition=npaodesc)
1840 22504 : ALLOCATE (qs_kind%pao_descriptors(npaodesc))
1841 11255 : DO ipaodesc = 1, npaodesc
1842 : CALL section_vals_val_get(pao_desc_section, keyword_name="BETA", i_rep_section=ipaodesc, &
1843 18 : 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 18 : 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 11255 : 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 11237 : keyword_name="ELEC_CONF", n_rep_val=i)
1853 11237 : IF (i > 0) THEN
1854 : CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
1855 4 : keyword_name="ELEC_CONF", i_vals=elec_conf)
1856 4 : 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 11237 : 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 11237 : keyword_name="POTENTIAL_FILE_NAME", c_val=potential_fn_kind)
1863 : CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
1864 11237 : keyword_name="POTENTIAL_TYPE", c_val=potential_type)
1865 : CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
1866 11237 : explicit=explicit, keyword_name="POTENTIAL", c_vals=tmpstringlist)
1867 11237 : IF (explicit) THEN
1868 10915 : IF (SIZE(tmpstringlist) == 1) THEN
1869 : ! old type of input: start of name defines type
1870 10807 : potential_name = tmpstringlist(1)
1871 10807 : IF (potential_type == "") THEN
1872 10807 : ipos = INDEX(potential_name, "-")
1873 10807 : IF (ipos > 1) THEN
1874 9505 : potential_type = potential_name(:ipos - 1)
1875 : ELSE
1876 1302 : potential_type = potential_name
1877 : END IF
1878 : END IF
1879 108 : ELSE IF (SIZE(tmpstringlist) == 2) THEN
1880 108 : potential_type = tmpstringlist(1)
1881 108 : potential_name = tmpstringlist(2)
1882 : ELSE
1883 0 : CPABORT("POTENTIAL input list is not correct")
1884 : END IF
1885 : END IF
1886 11237 : CALL uppercase(potential_type)
1887 :
1888 : ! Parse KG POTENTIAL
1889 : CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
1890 11237 : 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 11237 : 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 11237 : keyword_name="ECP_SEMI_LOCAL", l_val=ecp_semi_local)
1897 :
1898 : ! Assign atomic covalent radius
1899 11237 : qs_kind%covalent_radius = ptable(z)%covalent_radius*bohr
1900 : CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
1901 11237 : keyword_name="COVALENT_RADIUS", r_val=r)
1902 11237 : IF (r > 0.0_dp) qs_kind%covalent_radius = r
1903 :
1904 : ! Assign atomic van der Waals radius
1905 11237 : qs_kind%vdw_radius = ptable(z)%vdw_radius*bohr
1906 : CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
1907 11237 : keyword_name="VDW_RADIUS", r_val=r)
1908 11237 : 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 11237 : keyword_name="HARD_EXP_RADIUS")
1913 11237 : IF (i == 0) THEN
1914 11145 : IF (z == 1) THEN
1915 4908 : qs_kind%hard_radius = 1.2_dp
1916 : ELSE
1917 6237 : 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 92 : 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 11237 : keyword_name="RHO0_EXP_RADIUS")
1927 11237 : IF (i == 0) THEN
1928 11237 : qs_kind%hard0_radius = qs_kind%hard_radius
1929 : ELSE
1930 : CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
1931 0 : keyword_name="RHO0_EXP_RADIUS", r_val=qs_kind%hard0_radius)
1932 : END IF
1933 11237 : IF (qs_kind%hard_radius < qs_kind%hard0_radius) THEN
1934 0 : CPABORT("rc0 should be <= rc")
1935 : END IF
1936 :
1937 : CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
1938 11237 : keyword_name="MAX_RAD_LOCAL", r_val=qs_kind%max_rad_local)
1939 : CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
1940 11237 : keyword_name="LEBEDEV_GRID", i_val=qs_kind%ngrid_ang)
1941 11237 : IF (qs_kind%ngrid_ang <= 0) THEN
1942 0 : CPABORT("# point lebedev grid < 0")
1943 : END IF
1944 : CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
1945 11237 : keyword_name="RADIAL_GRID", i_val=qs_kind%ngrid_rad)
1946 11237 : IF (qs_kind%ngrid_rad <= 0) THEN
1947 0 : CPABORT("# point radial grid < 0")
1948 : END IF
1949 : CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
1950 11237 : keyword_name="GPW_TYPE", l_val=qs_kind%gpw_type_forced)
1951 : CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
1952 11237 : keyword_name="GHOST", l_val=qs_kind%ghost)
1953 : CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
1954 11237 : keyword_name="FLOATING_BASIS_CENTER", l_val=qs_kind%floating)
1955 : CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
1956 11237 : keyword_name="NO_OPTIMIZE", l_val=qs_kind%no_optimize)
1957 :
1958 : ! Magnetization
1959 : CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
1960 11237 : keyword_name="MAGNETIZATION", r_val=qs_kind%magnetization)
1961 : ! DFTB3 param
1962 : CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
1963 11237 : keyword_name="DFTB3_PARAM", r_val=qs_kind%dudq_dftb3)
1964 : CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
1965 11237 : keyword_name="LMAX_DFTB", i_val=qs_kind%lmax_dftb)
1966 :
1967 : ! MAOS
1968 : CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
1969 11237 : keyword_name="MAO", i_val=qs_kind%mao)
1970 :
1971 : ! Use monovalent pseudopotential
1972 : CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
1973 11237 : keyword_name="MONOVALENT", l_val=qs_kind%monovalent)
1974 11237 : IF (qs_kind%monovalent .AND. TRIM(potential_type) /= 'GTH') THEN
1975 0 : CPABORT("Monovalent pseudopotentials currently implemented for GTH only!")
1976 : END IF
1977 :
1978 : ! Read the BS subsection of the current atomic kind, if enabled
1979 11237 : NULLIFY (bs_section)
1980 : bs_section => section_vals_get_subs_vals(kind_section, "BS", &
1981 11237 : i_rep_section=k_rep)
1982 : section_enabled = .FALSE.
1983 : CALL section_vals_val_get(bs_section, "_SECTION_PARAMETERS_", &
1984 11237 : l_val=section_enabled)
1985 11237 : IF (section_enabled) THEN
1986 : ! test for conflict with magnetization
1987 66 : IF (qs_kind%magnetization /= 0.0_dp) THEN
1988 : CALL cp_abort(__LOCATION__, "BS Section is in conflict with non-zero magnetization "// &
1989 0 : "for this atom kind.")
1990 : END IF
1991 66 : qs_kind%bs_occupation = .TRUE.
1992 : !Alpha spin
1993 66 : NULLIFY (spin_section)
1994 66 : spin_section => section_vals_get_subs_vals(bs_section, "ALPHA")
1995 66 : CALL section_vals_get(spin_section, explicit=explicit)
1996 66 : IF (explicit) THEN
1997 66 : NULLIFY (add_el)
1998 : CALL section_vals_val_get(spin_section, &
1999 66 : keyword_name="NEL", i_vals=add_el)
2000 66 : CPASSERT(ASSOCIATED(add_el))
2001 198 : ALLOCATE (qs_kind%addel(SIZE(add_el), 2))
2002 358 : qs_kind%addel = 0
2003 146 : qs_kind%addel(1:SIZE(add_el), 1) = add_el(1:SIZE(add_el))
2004 66 : NULLIFY (add_el)
2005 : CALL section_vals_val_get(spin_section, &
2006 66 : keyword_name="L", i_vals=add_el)
2007 66 : CPASSERT(ASSOCIATED(add_el))
2008 66 : CPASSERT(SIZE(add_el) == SIZE(qs_kind%addel, 1))
2009 198 : ALLOCATE (qs_kind%laddel(SIZE(add_el), 2))
2010 358 : qs_kind%laddel = 0
2011 146 : qs_kind%laddel(1:SIZE(add_el), 1) = add_el(1:SIZE(add_el))
2012 198 : ALLOCATE (qs_kind%naddel(SIZE(add_el), 2))
2013 358 : qs_kind%naddel = 0
2014 66 : NULLIFY (add_el)
2015 : CALL section_vals_val_get(spin_section, &
2016 66 : keyword_name="N", n_rep_val=i)
2017 66 : IF (i > 0) THEN
2018 : CALL section_vals_val_get(spin_section, &
2019 66 : keyword_name="N", i_vals=add_el)
2020 66 : IF (SIZE(add_el) == SIZE(qs_kind%addel, 1)) THEN
2021 146 : qs_kind%naddel(1:SIZE(add_el), 1) = add_el(1:SIZE(add_el))
2022 : END IF
2023 : END IF
2024 : END IF
2025 : ! Beta spin
2026 66 : NULLIFY (spin_section)
2027 66 : spin_section => section_vals_get_subs_vals(bs_section, "BETA")
2028 66 : CALL section_vals_get(spin_section, explicit=explicit)
2029 66 : IF (explicit) THEN
2030 66 : NULLIFY (add_el)
2031 : CALL section_vals_val_get(spin_section, &
2032 66 : keyword_name="NEL", i_vals=add_el)
2033 66 : CPASSERT(SIZE(add_el) == SIZE(qs_kind%addel, 1))
2034 146 : qs_kind%addel(1:SIZE(add_el), 2) = add_el(1:SIZE(add_el))
2035 358 : qs_kind%addel(:, :) = qs_kind%addel(:, :)
2036 66 : NULLIFY (add_el)
2037 : CALL section_vals_val_get(spin_section, &
2038 66 : keyword_name="L", i_vals=add_el)
2039 66 : CPASSERT(SIZE(add_el) == SIZE(qs_kind%addel, 1))
2040 146 : qs_kind%laddel(1:SIZE(add_el), 2) = add_el(1:SIZE(add_el))
2041 :
2042 : CALL section_vals_val_get(spin_section, &
2043 66 : keyword_name="N", n_rep_val=i)
2044 66 : IF (i > 0) THEN
2045 66 : NULLIFY (add_el)
2046 : CALL section_vals_val_get(spin_section, &
2047 66 : keyword_name="N", i_vals=add_el)
2048 66 : IF (SIZE(add_el) == SIZE(qs_kind%addel, 1)) THEN
2049 146 : qs_kind%naddel(1:SIZE(add_el), 2) = add_el(1:SIZE(add_el))
2050 : END IF
2051 : END IF
2052 : END IF
2053 : END IF
2054 :
2055 : ! Read the DFT+U subsection of the current atomic kind, if enabled
2056 :
2057 11237 : NULLIFY (dft_plus_u_section)
2058 : dft_plus_u_section => section_vals_get_subs_vals(kind_section, &
2059 : subsection_name="DFT_PLUS_U", &
2060 11237 : i_rep_section=k_rep)
2061 : section_enabled = .FALSE.
2062 : CALL section_vals_val_get(dft_plus_u_section, &
2063 : keyword_name="_SECTION_PARAMETERS_", &
2064 11237 : l_val=section_enabled)
2065 146081 : IF (section_enabled) THEN
2066 36 : ALLOCATE (qs_kind%dft_plus_u)
2067 36 : NULLIFY (qs_kind%dft_plus_u%nelec)
2068 36 : NULLIFY (qs_kind%dft_plus_u%orbitals)
2069 : CALL section_vals_val_get(dft_plus_u_section, &
2070 : keyword_name="L", &
2071 36 : i_val=l)
2072 36 : qs_kind%dft_plus_u%l = l
2073 : #if defined(__SIRIUS)
2074 : CALL section_vals_val_get(dft_plus_u_section, &
2075 : keyword_name="N", &
2076 36 : i_val=nu)
2077 36 : qs_kind%dft_plus_u%n = nu
2078 :
2079 : CALL section_vals_val_get(dft_plus_u_section, &
2080 : keyword_name="U", &
2081 : r_val=qs_kind%dft_plus_u%U, &
2082 36 : explicit=explicit_U)
2083 :
2084 : CALL section_vals_val_get(dft_plus_u_section, &
2085 : keyword_name="J", &
2086 : r_val=qs_kind%dft_plus_u%J, &
2087 36 : explicit=explicit_J)
2088 :
2089 : CALL section_vals_val_get(dft_plus_u_section, &
2090 : keyword_name="alpha", &
2091 36 : r_val=qs_kind%dft_plus_u%alpha)
2092 :
2093 : CALL section_vals_val_get(dft_plus_u_section, &
2094 : keyword_name="beta", &
2095 36 : r_val=qs_kind%dft_plus_u%beta)
2096 :
2097 : CALL section_vals_val_get(dft_plus_u_section, &
2098 : keyword_name="J0", &
2099 36 : r_val=qs_kind%dft_plus_u%J0)
2100 :
2101 : CALL section_vals_val_get(dft_plus_u_section, &
2102 : keyword_name="occupation", &
2103 36 : r_val=qs_kind%dft_plus_u%occupation)
2104 : #else
2105 : nu = 0
2106 : #endif
2107 :
2108 : CALL section_vals_val_get(dft_plus_u_section, &
2109 : keyword_name="U_MINUS_J", &
2110 : r_val=qs_kind%dft_plus_u%u_minus_j_target, &
2111 36 : explicit=explicit_u_m_j)
2112 :
2113 36 : IF ((explicit_U .OR. explicit_J) .AND. explicit_u_m_j) THEN
2114 0 : CPABORT("DFT+U| specifying U or J and U_MINUS_J parameters are mutually exclusive.")
2115 : END IF
2116 :
2117 : CALL section_vals_val_get(dft_plus_u_section, &
2118 : keyword_name="U_RAMPING", &
2119 36 : r_val=qs_kind%dft_plus_u%u_ramping)
2120 : CALL section_vals_val_get(dft_plus_u_section, &
2121 : keyword_name="INIT_U_RAMPING_EACH_SCF", &
2122 36 : l_val=qs_kind%dft_plus_u%init_u_ramping_each_scf)
2123 36 : IF (qs_kind%dft_plus_u%u_ramping > 0.0_dp) THEN
2124 8 : qs_kind%dft_plus_u%u_minus_j = 0.0_dp
2125 : ELSE
2126 28 : qs_kind%dft_plus_u%u_minus_j = qs_kind%dft_plus_u%u_minus_j_target
2127 : END IF
2128 : CALL section_vals_val_get(dft_plus_u_section, &
2129 : keyword_name="EPS_U_RAMPING", &
2130 36 : r_val=qs_kind%dft_plus_u%eps_u_ramping)
2131 :
2132 36 : NULLIFY (enforce_occupation_section)
2133 : enforce_occupation_section => section_vals_get_subs_vals(dft_plus_u_section, &
2134 36 : subsection_name="ENFORCE_OCCUPATION")
2135 : subsection_enabled = .FALSE.
2136 : CALL section_vals_val_get(enforce_occupation_section, &
2137 : keyword_name="_SECTION_PARAMETERS_", &
2138 36 : l_val=subsection_enabled)
2139 36 : IF (subsection_enabled) THEN
2140 4 : NULLIFY (nelec)
2141 : CALL section_vals_val_get(enforce_occupation_section, &
2142 : keyword_name="NELEC", &
2143 4 : r_vals=nelec)
2144 4 : nspin = SIZE(nelec)
2145 12 : ALLOCATE (qs_kind%dft_plus_u%nelec(nspin))
2146 8 : qs_kind%dft_plus_u%nelec(:) = nelec(:)
2147 4 : NULLIFY (orbitals)
2148 : CALL section_vals_val_get(enforce_occupation_section, &
2149 : keyword_name="ORBITALS", &
2150 4 : i_vals=orbitals)
2151 4 : norbitals = SIZE(orbitals)
2152 4 : IF (norbitals <= 0 .OR. norbitals > 2*l + 1) THEN
2153 : CALL cp_abort(__LOCATION__, "DFT+U| Invalid number of ORBITALS specified: "// &
2154 0 : "1 to 2*L+1 integer numbers are expected")
2155 : END IF
2156 12 : ALLOCATE (qs_kind%dft_plus_u%orbitals(norbitals))
2157 16 : qs_kind%dft_plus_u%orbitals(:) = orbitals(:)
2158 4 : NULLIFY (orbitals)
2159 16 : DO m = 1, norbitals
2160 12 : IF (qs_kind%dft_plus_u%orbitals(m) > l) THEN
2161 0 : CPABORT("DFT+U| Invalid orbital magnetic quantum number specified: m > l")
2162 : END IF
2163 12 : IF (qs_kind%dft_plus_u%orbitals(m) < -l) THEN
2164 0 : CPABORT("DFT+U| Invalid orbital magnetic quantum number specified: m < -l")
2165 : END IF
2166 52 : DO j = 1, norbitals
2167 48 : IF (j /= m) THEN
2168 24 : IF (qs_kind%dft_plus_u%orbitals(j) == qs_kind%dft_plus_u%orbitals(m)) THEN
2169 0 : CPABORT("DFT+U| An orbital magnetic quantum number was specified twice")
2170 : END IF
2171 : END IF
2172 : END DO
2173 : END DO
2174 : CALL section_vals_val_get(enforce_occupation_section, &
2175 : keyword_name="EPS_SCF", &
2176 4 : r_val=qs_kind%dft_plus_u%eps_scf)
2177 : CALL section_vals_val_get(enforce_occupation_section, &
2178 : keyword_name="MAX_SCF", &
2179 4 : i_val=i)
2180 4 : qs_kind%dft_plus_u%max_scf = MAX(-1, i)
2181 : CALL section_vals_val_get(enforce_occupation_section, &
2182 : keyword_name="SMEAR", &
2183 4 : l_val=qs_kind%dft_plus_u%smear)
2184 : END IF ! subsection enabled
2185 : END IF ! section enabled
2186 :
2187 : END IF
2188 :
2189 : ! Allocate and initialise the orbital basis set data set structure
2190 16447 : CALL init_orbital_pointers(5) ! debug the SUN optimizer
2191 :
2192 : ! BASIS and POTENTIAL read only when strictly necessary otherwise, even if not used
2193 : ! we just print misleading informations
2194 16447 : explicit_basis = .FALSE.
2195 16447 : IF (k_rep > 0) THEN
2196 : basis_section => section_vals_get_subs_vals(kind_section, "BASIS", i_rep_section=k_rep, &
2197 11237 : can_return_null=.TRUE.)
2198 11237 : CALL section_vals_get(basis_section, explicit=explicit_basis)
2199 : END IF
2200 :
2201 16447 : explicit_potential = .FALSE.
2202 16447 : IF (k_rep > 0) THEN
2203 : potential_section => section_vals_get_subs_vals(kind_section, "POTENTIAL", &
2204 11237 : i_rep_section=k_rep, can_return_null=.TRUE.)
2205 11237 : CALL section_vals_get(potential_section, explicit=explicit_potential)
2206 : END IF
2207 :
2208 16447 : explicit_kgpot = .FALSE.
2209 16447 : IF (k_rep > 0) THEN
2210 : kgpot_section => section_vals_get_subs_vals(kind_section, "KG_POTENTIAL", &
2211 11237 : i_rep_section=k_rep, can_return_null=.TRUE.)
2212 11237 : CALL section_vals_get(kgpot_section, explicit=explicit_kgpot)
2213 : END IF
2214 :
2215 18691 : SELECT CASE (method_id)
2216 : CASE (do_method_rm1, do_method_am1, do_method_mndo, do_method_pdg, do_method_pm3, do_method_pm6, &
2217 : do_method_pm6fm, do_method_mndod, do_method_pnnl)
2218 : ! Allocate all_potential
2219 2244 : CALL allocate_potential(qs_kind%all_potential)
2220 2244 : CALL set_default_all_potential(qs_kind%all_potential, z, zeff_correction)
2221 2244 : CALL get_qs_kind(qs_kind, elec_conf=elec_conf)
2222 2244 : IF (.NOT. ASSOCIATED(elec_conf)) THEN
2223 2244 : CALL get_potential(potential=qs_kind%all_potential, elec_conf=elec_conf)
2224 2244 : CALL set_qs_kind(qs_kind, elec_conf=elec_conf)
2225 : END IF
2226 2244 : CPASSERT(.NOT. qs_kind%floating)
2227 2244 : IF (qs_kind%ghost) THEN
2228 0 : CALL get_qs_kind(qs_kind=qs_kind, elec_conf=elec_conf)
2229 0 : elec_conf(:) = 0
2230 : CALL get_potential(potential=qs_kind%all_potential, &
2231 0 : elec_conf=elec_conf)
2232 0 : elec_conf(:) = 0
2233 : CALL set_potential(potential=qs_kind%all_potential, &
2234 : zeff=0.0_dp, &
2235 0 : zeff_correction=0.0_dp)
2236 : END IF
2237 :
2238 : ! Basis set (Parameters)
2239 : ! Setup proper semiempirical parameters
2240 2244 : check = .NOT. ASSOCIATED(qs_kind%se_parameter)
2241 2244 : CPASSERT(check)
2242 2244 : CALL semi_empirical_create(qs_kind%se_parameter)
2243 : ! Check if we allow p-orbitals on H
2244 440 : SELECT CASE (z)
2245 : CASE (1)
2246 2244 : IF (k_rep > 0) THEN
2247 : CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
2248 52 : keyword_name="SE_P_ORBITALS_ON_H", l_val=qs_kind%se_parameter%p_orbitals_on_h)
2249 : END IF
2250 : CASE DEFAULT
2251 : ! No special cases for other elements..
2252 : END SELECT
2253 : ! Set default parameters
2254 2244 : CALL section_vals_val_get(dft_section, "QS%SE%STO_NG", i_val=ngauss)
2255 2244 : CALL se_param_set_default(qs_kind%se_parameter, z, method_id)
2256 2244 : NULLIFY (tmp_basis_set)
2257 2244 : CALL init_se_param(qs_kind%se_parameter, tmp_basis_set, ngauss)
2258 2244 : CALL add_basis_set_to_container(qs_kind%basis_sets, tmp_basis_set, "ORB")
2259 : CALL init_potential(qs_kind%all_potential, itype="BARE", &
2260 2244 : zeff=qs_kind%se_parameter%zeff, zeff_correction=zeff_correction)
2261 2244 : qs_kind%se_parameter%zeff = qs_kind%se_parameter%zeff - zeff_correction
2262 :
2263 2244 : check = ((potential_name /= '') .OR. explicit_potential) .AND. .NOT. silent
2264 : IF (check) THEN
2265 : CALL cp_warn(__LOCATION__, &
2266 : "Information provided in the input file regarding POTENTIAL for KIND <"// &
2267 80 : TRIM(qs_kind%name)//"> will be ignored!")
2268 : END IF
2269 :
2270 2244 : check = ((k_rep > 0) .OR. explicit_basis) .AND. .NOT. silent
2271 618 : IF (check) THEN
2272 : CALL cp_warn(__LOCATION__, &
2273 : "Information provided in the input file regarding BASIS for KIND <"// &
2274 116 : TRIM(qs_kind%name)//"> will be ignored!")
2275 : END IF
2276 :
2277 : CASE (do_method_dftb)
2278 : ! Allocate all_potential
2279 618 : CALL allocate_potential(qs_kind%all_potential)
2280 618 : CALL set_default_all_potential(qs_kind%all_potential, z, zeff_correction)
2281 618 : CALL get_qs_kind(qs_kind, elec_conf=elec_conf)
2282 618 : IF (.NOT. ASSOCIATED(elec_conf)) THEN
2283 618 : CALL get_potential(potential=qs_kind%all_potential, elec_conf=elec_conf)
2284 618 : CALL set_qs_kind(qs_kind, elec_conf=elec_conf)
2285 : END IF
2286 618 : CPASSERT(.NOT. qs_kind%floating)
2287 618 : IF (qs_kind%ghost) THEN
2288 0 : CALL get_qs_kind(qs_kind=qs_kind, elec_conf=elec_conf)
2289 0 : elec_conf(:) = 0
2290 : CALL get_potential(potential=qs_kind%all_potential, &
2291 0 : elec_conf=elec_conf)
2292 0 : elec_conf(:) = 0
2293 : CALL set_potential(potential=qs_kind%all_potential, &
2294 : zeff=0.0_dp, &
2295 0 : zeff_correction=0.0_dp)
2296 : END IF
2297 :
2298 618 : check = ((potential_name /= '') .OR. explicit_potential) .AND. .NOT. silent
2299 : IF (check) THEN
2300 : CALL cp_warn(__LOCATION__, &
2301 : "Information provided in the input file regarding POTENTIAL for KIND <"// &
2302 0 : TRIM(qs_kind%name)//"> will be ignored!")
2303 : END IF
2304 :
2305 618 : check = ((k_rep > 0) .OR. explicit_basis) .AND. .NOT. silent
2306 2562 : IF (check) THEN
2307 : CALL cp_warn(__LOCATION__, &
2308 : "Information provided in the input file regarding BASIS for KIND <"// &
2309 90 : TRIM(qs_kind%name)//"> will be ignored!")
2310 : END IF
2311 :
2312 : CASE (do_method_xtb)
2313 : ! Allocate all_potential
2314 2562 : CALL allocate_potential(qs_kind%all_potential)
2315 2562 : CALL set_default_all_potential(qs_kind%all_potential, z, zeff_correction)
2316 2562 : CALL get_qs_kind(qs_kind, elec_conf=elec_conf)
2317 2562 : IF (.NOT. ASSOCIATED(elec_conf)) THEN
2318 2562 : CALL get_potential(potential=qs_kind%all_potential, elec_conf=elec_conf)
2319 2562 : CALL set_qs_kind(qs_kind, elec_conf=elec_conf)
2320 : END IF
2321 2562 : CPASSERT(.NOT. qs_kind%floating)
2322 2562 : IF (qs_kind%ghost) THEN
2323 0 : CALL get_qs_kind(qs_kind=qs_kind, elec_conf=elec_conf)
2324 0 : elec_conf(:) = 0
2325 : CALL get_potential(potential=qs_kind%all_potential, &
2326 0 : elec_conf=elec_conf)
2327 0 : elec_conf(:) = 0
2328 : CALL set_potential(potential=qs_kind%all_potential, &
2329 : zeff=0.0_dp, &
2330 0 : zeff_correction=0.0_dp)
2331 : END IF
2332 :
2333 2562 : check = ((potential_name /= '') .OR. explicit_potential) .AND. .NOT. silent
2334 : IF (check) THEN
2335 : CALL cp_warn(__LOCATION__, &
2336 : "Information provided in the input file regarding POTENTIAL for KIND <"// &
2337 0 : TRIM(qs_kind%name)//"> will be ignored!")
2338 : END IF
2339 :
2340 2562 : check = ((k_rep > 0) .OR. explicit_basis) .AND. .NOT. silent
2341 26 : IF (check) THEN
2342 : CALL cp_warn(__LOCATION__, &
2343 : "Information provided in the input file regarding BASIS for KIND <"// &
2344 0 : TRIM(qs_kind%name)//"> will be ignored!")
2345 : END IF
2346 :
2347 : CASE (do_method_pw)
2348 : ! PW DFT
2349 : ! Allocate and initialise the potential data set structure
2350 26 : IF (potential_name /= '') THEN
2351 26 : SELECT CASE (TRIM(potential_type))
2352 : CASE ("ALL", "ECP", "CNEO")
2353 : CALL cp_abort(__LOCATION__, &
2354 : "PW DFT calculations only with potential type UPF or GTH possible."// &
2355 : " <"//TRIM(potential_type)//"> was specified "// &
2356 0 : "for the atomic kind <"//TRIM(qs_kind%name))
2357 : CASE ("GTH")
2358 6 : IF (potential_fn_kind == "-") THEN
2359 6 : CALL section_vals_val_get(dft_section, "POTENTIAL_FILE_NAME", c_val=potential_file_name)
2360 : ELSE
2361 0 : potential_file_name = potential_fn_kind
2362 : END IF
2363 6 : CALL allocate_potential(qs_kind%gth_potential)
2364 : CALL read_potential(qs_kind%element_symbol, potential_name, &
2365 : qs_kind%gth_potential, zeff_correction, para_env, &
2366 : potential_file_name, potential_section, update_input, &
2367 6 : monovalent=qs_kind%monovalent)
2368 6 : CALL set_potential(qs_kind%gth_potential, z=z)
2369 6 : CALL get_qs_kind(qs_kind, elec_conf=elec_conf)
2370 6 : IF (.NOT. ASSOCIATED(elec_conf)) THEN
2371 6 : CALL get_potential(potential=qs_kind%gth_potential, elec_conf=elec_conf)
2372 6 : CALL set_qs_kind(qs_kind, elec_conf=elec_conf)
2373 : ELSE
2374 0 : CALL set_potential(potential=qs_kind%gth_potential, elec_conf=elec_conf)
2375 : END IF
2376 : CASE ("UPF")
2377 2100 : ALLOCATE (qs_kind%upf_potential)
2378 20 : qs_kind%upf_potential%zion = 0
2379 20 : qs_kind%upf_potential%filename = ADJUSTL(TRIM(potential_name))
2380 20 : CALL get_qs_kind(qs_kind, elec_conf=elec_conf)
2381 20 : IF (.NOT. ASSOCIATED(elec_conf)) THEN
2382 20 : CALL set_qs_kind(qs_kind, elec_conf=qs_kind%upf_potential%econf)
2383 : END IF
2384 : CASE DEFAULT
2385 : CALL cp_abort(__LOCATION__, &
2386 : "An invalid potential type <"// &
2387 : TRIM(potential_type)//"> was specified "// &
2388 : "for the atomic kind <"// &
2389 26 : TRIM(qs_kind%name))
2390 : END SELECT
2391 : ELSE
2392 : CALL cp_abort(__LOCATION__, &
2393 : "No potential type was defined for the "// &
2394 0 : "atomic kind <"//TRIM(qs_kind%name)//">")
2395 : END IF
2396 :
2397 : CASE DEFAULT
2398 :
2399 : ! set ngauss for STO expansion
2400 10997 : CALL section_vals_val_get(dft_section, "QS%STO_NG", i_val=ngauss)
2401 : ! Allocate and initialise the basis set data set structure
2402 : ! first external basis sets
2403 24117 : DO i = 1, nb_rep
2404 26236 : SELECT CASE (basis_set_form(i))
2405 : CASE ("GTO")
2406 13116 : NULLIFY (tmp_basis_set)
2407 13116 : CALL allocate_gto_basis_set(tmp_basis_set)
2408 : CALL read_gto_basis_set(qs_kind%element_symbol, basis_set_name(i), &
2409 13116 : tmp_basis_set, para_env, dft_section)
2410 : CASE ("STO")
2411 4 : NULLIFY (sto_basis_set)
2412 4 : CALL allocate_sto_basis_set(sto_basis_set)
2413 : CALL read_sto_basis_set(qs_kind%element_symbol, basis_set_name(i), &
2414 4 : sto_basis_set, para_env, dft_section)
2415 4 : NULLIFY (tmp_basis_set)
2416 4 : CALL create_gto_from_sto_basis(sto_basis_set, tmp_basis_set, ngauss)
2417 4 : CALL deallocate_sto_basis_set(sto_basis_set)
2418 : CASE DEFAULT
2419 : CALL cp_abort(__LOCATION__, &
2420 : "Invalid basis set form "//TRIM(basis_set_form(i))// &
2421 13120 : "for atomic kind <"//TRIM(qs_kind%name)//">")
2422 : END SELECT
2423 13120 : tmp = basis_set_type(i)
2424 13120 : CALL uppercase(tmp)
2425 24117 : CALL add_basis_set_to_container(qs_kind%basis_sets, tmp_basis_set, tmp)
2426 : END DO
2427 : ! now explicit basis sets
2428 10997 : IF (explicit_basis) THEN
2429 162 : CALL section_vals_get(basis_section, n_repetition=nexp)
2430 324 : DO i = 1, nexp
2431 162 : NULLIFY (tmp_basis_set)
2432 162 : CALL allocate_gto_basis_set(tmp_basis_set)
2433 : CALL read_gto_basis_set(qs_kind%element_symbol, basis_type, &
2434 162 : tmp_basis_set, basis_section, i, dft_section)
2435 162 : tmp = basis_type
2436 162 : CALL uppercase(tmp)
2437 324 : CALL add_basis_set_to_container(qs_kind%basis_sets, tmp_basis_set, tmp)
2438 : END DO
2439 : END IF
2440 : ! combine multiple basis sets
2441 230937 : DO i = 1, SIZE(qs_kind%basis_sets)
2442 219940 : NULLIFY (tmp_basis_set)
2443 : CALL get_basis_from_container(qs_kind%basis_sets, basis_set=tmp_basis_set, &
2444 219940 : inumbas=i, basis_type=basis_type)
2445 219940 : IF (basis_type == "") CYCLE
2446 13282 : jj = i
2447 263225 : DO j = i + 1, SIZE(qs_kind%basis_sets)
2448 249943 : jj = jj + 1
2449 249943 : NULLIFY (sup_basis_set)
2450 : CALL get_basis_from_container(qs_kind%basis_sets, basis_set=sup_basis_set, &
2451 249943 : inumbas=jj, basis_type=tmp)
2452 263225 : IF (basis_type == tmp) THEN
2453 : ! we found a match, combine the basis sets and delete the second
2454 0 : CALL combine_basis_sets(tmp_basis_set, sup_basis_set)
2455 0 : CALL remove_basis_from_container(qs_kind%basis_sets, jj)
2456 0 : jj = jj - 1
2457 : END IF
2458 : END DO
2459 230937 : NULLIFY (sup_basis_set)
2460 : END DO
2461 :
2462 : ! check that we have an orbital basis set
2463 10997 : nobasis = .TRUE.
2464 230937 : DO i = 1, SIZE(qs_kind%basis_sets)
2465 219940 : NULLIFY (tmp_basis_set)
2466 : CALL get_basis_from_container(qs_kind%basis_sets, basis_set=tmp_basis_set, &
2467 219940 : inumbas=i, basis_type=basis_type)
2468 230937 : IF (basis_type == "ORB") nobasis = .FALSE.
2469 : END DO
2470 10997 : IF (nobasis) THEN
2471 : CALL cp_abort(__LOCATION__, &
2472 : "No basis set type was defined for the "// &
2473 0 : "atomic kind <"//TRIM(qs_kind%name)//">")
2474 : END IF
2475 :
2476 : ! If Ghost atom we don't need to allocate/initialize anything connected to POTENTIAL
2477 10997 : IF (qs_kind%ghost .OR. qs_kind%floating) THEN
2478 150 : IF (ASSOCIATED(qs_kind%elec_conf)) qs_kind%elec_conf = 0
2479 : ELSE
2480 : ! Allocate and initialise the potential data set structure
2481 10847 : IF ((potential_name /= '') .OR. explicit_potential) THEN
2482 : ! determine the pseudopotential file to search
2483 10847 : IF (potential_fn_kind == "-") THEN
2484 10837 : CALL section_vals_val_get(dft_section, "POTENTIAL_FILE_NAME", c_val=potential_file_name)
2485 : ELSE
2486 10 : potential_file_name = potential_fn_kind
2487 : END IF
2488 : !
2489 12137 : SELECT CASE (TRIM(potential_type))
2490 : CASE ("ALL")
2491 1290 : CALL allocate_potential(qs_kind%all_potential)
2492 : CALL read_potential(qs_kind%element_symbol, potential_name, &
2493 : qs_kind%all_potential, zeff_correction, para_env, &
2494 1290 : potential_file_name, potential_section, update_input)
2495 1290 : CALL set_potential(qs_kind%all_potential, z=z)
2496 1290 : CALL get_qs_kind(qs_kind, elec_conf=elec_conf)
2497 1290 : IF (.NOT. ASSOCIATED(elec_conf)) THEN
2498 1290 : CALL get_potential(potential=qs_kind%all_potential, elec_conf=elec_conf)
2499 1290 : CALL set_qs_kind(qs_kind, elec_conf=elec_conf)
2500 : ELSE
2501 0 : CALL set_potential(potential=qs_kind%all_potential, elec_conf=elec_conf)
2502 : END IF
2503 : CASE ("GTH")
2504 9461 : CALL allocate_potential(qs_kind%gth_potential)
2505 : CALL read_potential(qs_kind%element_symbol, potential_name, &
2506 : qs_kind%gth_potential, zeff_correction, para_env, &
2507 : potential_file_name, potential_section, update_input, &
2508 9461 : monovalent=qs_kind%monovalent)
2509 9461 : CALL set_potential(qs_kind%gth_potential, z=z)
2510 9461 : CALL get_qs_kind(qs_kind, elec_conf=elec_conf)
2511 9461 : IF (.NOT. ASSOCIATED(elec_conf)) THEN
2512 9457 : CALL get_potential(potential=qs_kind%gth_potential, elec_conf=elec_conf)
2513 9457 : CALL set_qs_kind(qs_kind, elec_conf=elec_conf)
2514 : ELSE
2515 4 : CALL set_potential(potential=qs_kind%gth_potential, elec_conf=elec_conf)
2516 : END IF
2517 : CASE ("ECP")
2518 76 : CALL allocate_potential(qs_kind%sgp_potential)
2519 76 : CALL get_potential(qs_kind%sgp_potential, description=description)
2520 : CALL read_ecp_potential(ptable(z)%symbol, ecppot, &
2521 76 : potential_name, potential_file_name, potential_section)
2522 76 : IF (ecp_semi_local) THEN
2523 76 : description(1) = "Semi-local Gaussian pseudopotential "
2524 76 : description(2) = "ECP "//TRIM(potential_name)
2525 76 : description(3) = "LIBGRPP: A. V. Oleynichenko et al., Symmetry 15 197 2023"
2526 : description(4) = " "
2527 : ELSE
2528 0 : description(4) = "ECP "//TRIM(potential_name)
2529 : END IF
2530 : CALL set_potential(qs_kind%sgp_potential, name=ecppot%pname, description=description, &
2531 : zeff=ecppot%zion, z=z, ecp_local=.TRUE., ecp_semi_local=ecp_semi_local, &
2532 : nloc=ecppot%nloc, nrloc=ecppot%nrloc, aloc=ecppot%aloc, bloc=ecppot%bloc, &
2533 76 : has_nlcc=.FALSE.)
2534 : CALL set_potential(qs_kind%sgp_potential, sl_lmax=ecppot%lmax, &
2535 76 : npot=ecppot%npot, nrpot=ecppot%nrpot, apot=ecppot%apot, bpot=ecppot%bpot)
2536 : ! convert PP
2537 76 : IF (.NOT. ecp_semi_local) THEN
2538 0 : CPABORT("ECPs are only well tested in their semi-local form")
2539 0 : CALL get_qs_kind(qs_kind, basis_set=orb_basis_set)
2540 0 : CALL sgp_construction(sgp_pot=sgppot, ecp_pot=ecppot, orb_basis=orb_basis_set, error=error)
2541 0 : IF (iounit > 0 .AND. .NOT. silent) THEN
2542 0 : WRITE (iounit, "(/,T2,'PP Transformation for ',A)") TRIM(ecppot%pname)
2543 0 : IF (sgppot%has_local) THEN
2544 0 : WRITE (iounit, "(T8,'Accuracy for local part:',T41,F10.3,'%',T61,F20.12)") error(4), error(1)
2545 : END IF
2546 0 : IF (sgppot%has_nonlocal) THEN
2547 0 : WRITE (iounit, "(T8,'Accuracy for nonlocal part:',T41,F10.3,'%',T61,F20.12)") error(5), error(2)
2548 : END IF
2549 0 : IF (sgppot%has_nlcc) THEN
2550 0 : WRITE (iounit, "(T8,'Accuracy for NLCC density:',T61,F20.12)") error(3)
2551 : END IF
2552 : END IF
2553 : END IF
2554 76 : IF (sgppot%has_nonlocal) THEN
2555 : CALL set_potential(qs_kind%sgp_potential, n_nonlocal=sgppot%n_nonlocal, lmax=sgppot%lmax, &
2556 0 : is_nonlocal=sgppot%is_nonlocal)
2557 0 : nnl = sgppot%n_nonlocal
2558 0 : nppnl = 0
2559 0 : DO l = 0, sgppot%lmax
2560 0 : nppnl = nppnl + nnl*nco(l)
2561 : END DO
2562 0 : l = sgppot%lmax
2563 0 : ALLOCATE (a_nl(nnl), h_nl(nnl, 0:l), c_nl(nnl, nnl, 0:l))
2564 0 : a_nl(:) = sgppot%a_nonlocal(:)
2565 0 : h_nl(:, :) = sgppot%h_nonlocal(:, :)
2566 0 : DO l = 0, sgppot%lmax
2567 0 : c_nl(:, :, l) = sgppot%c_nonlocal(:, :, l)*SQRT(2._dp*l + 1.0_dp)
2568 : END DO
2569 0 : CALL set_potential(qs_kind%sgp_potential, nppnl=nppnl, a_nonlocal=a_nl, h_nonlocal=h_nl, c_nonlocal=c_nl)
2570 : ELSE
2571 76 : CALL set_potential(qs_kind%sgp_potential, n_nonlocal=0, lmax=-1, is_nonlocal=sgppot%is_nonlocal)
2572 76 : CALL set_potential(qs_kind%sgp_potential, nppnl=0)
2573 : END IF
2574 : !
2575 76 : CPASSERT(.NOT. sgppot%has_local)
2576 76 : CPASSERT(.NOT. sgppot%has_nlcc)
2577 : ! core
2578 76 : rc = 0.5_dp*qs_kind%covalent_radius*angstrom
2579 76 : rc = MAX(rc, 0.2_dp)
2580 76 : rc = MIN(rc, 1.0_dp)
2581 76 : alpha = 1.0_dp/(2.0_dp*rc**2)
2582 76 : ccore = ecppot%zion*SQRT((alpha/pi)**3)
2583 : CALL set_potential(qs_kind%sgp_potential, alpha_core_charge=alpha, ccore_charge=ccore, &
2584 76 : core_charge_radius=rc)
2585 76 : CALL atom_sgp_release(sgppot)
2586 76 : CALL get_qs_kind(qs_kind, elec_conf=elec_conf)
2587 76 : IF (.NOT. ASSOCIATED(elec_conf)) THEN
2588 76 : CALL set_qs_kind(qs_kind, elec_conf=ecppot%econf)
2589 : END IF
2590 76 : CALL get_qs_kind(qs_kind, elec_conf=elec_conf)
2591 76 : CALL set_potential(qs_kind%sgp_potential, elec_conf=elec_conf)
2592 : CASE ("UPF")
2593 12 : CALL allocate_potential(qs_kind%sgp_potential)
2594 12 : CALL get_potential(qs_kind%sgp_potential, description=description)
2595 12 : description(4) = "UPF "//TRIM(potential_name)
2596 12 : CALL atom_read_upf(upfpot, potential_name)
2597 : CALL set_potential(qs_kind%sgp_potential, name=upfpot%pname, description=description, &
2598 12 : zeff=upfpot%zion, z=z, has_nlcc=upfpot%core_correction)
2599 : ! convert pp
2600 12 : CALL sgp_construction(sgp_pot=sgppot, upf_pot=upfpot, error=error)
2601 12 : IF (iounit > 0 .AND. .NOT. silent) THEN
2602 6 : WRITE (iounit, "(/,T2,'PP Transformation for ',A)") TRIM(upfpot%pname)
2603 6 : IF (sgppot%has_local) THEN
2604 6 : WRITE (iounit, "(T8,'Accuracy for local part:',T61,F20.12)") error(1)
2605 : END IF
2606 6 : IF (sgppot%has_nonlocal) THEN
2607 3 : WRITE (iounit, "(T8,'Accuracy for nonlocal part:',T61,F20.12)") error(2)
2608 : END IF
2609 6 : IF (sgppot%has_nlcc) THEN
2610 0 : WRITE (iounit, "(T8,'Accuracy for NLCC density:',T61,F20.12)") error(3)
2611 : END IF
2612 : END IF
2613 12 : IF (sgppot%has_nonlocal) THEN
2614 : CALL set_potential(qs_kind%sgp_potential, n_nonlocal=sgppot%n_nonlocal, lmax=sgppot%lmax, &
2615 6 : is_nonlocal=sgppot%is_nonlocal)
2616 6 : nnl = sgppot%n_nonlocal
2617 6 : nppnl = 0
2618 12 : DO l = 0, sgppot%lmax
2619 12 : nppnl = nppnl + nnl*nco(l)
2620 : END DO
2621 6 : l = sgppot%lmax
2622 60 : ALLOCATE (a_nl(nnl), h_nl(nnl, 0:l), c_nl(nnl, nnl, 0:l))
2623 54 : a_nl(:) = sgppot%a_nonlocal(:)
2624 60 : h_nl(:, :) = sgppot%h_nonlocal(:, :)
2625 444 : c_nl(:, :, :) = sgppot%c_nonlocal(:, :, :)
2626 6 : CALL set_potential(qs_kind%sgp_potential, nppnl=nppnl, a_nonlocal=a_nl, h_nonlocal=h_nl, c_nonlocal=c_nl)
2627 : ELSE
2628 6 : CALL set_potential(qs_kind%sgp_potential, n_nonlocal=0, lmax=-1, is_nonlocal=sgppot%is_nonlocal)
2629 6 : CALL set_potential(qs_kind%sgp_potential, nppnl=0)
2630 : END IF
2631 12 : CPASSERT(sgppot%has_local)
2632 : ! core
2633 12 : rc = sgppot%ac_local
2634 12 : alpha = 1.0_dp/(2.0_dp*rc**2)
2635 12 : ccore = upfpot%zion*SQRT((alpha/pi)**3)
2636 : CALL set_potential(qs_kind%sgp_potential, alpha_core_charge=alpha, ccore_charge=ccore, &
2637 12 : core_charge_radius=rc)
2638 : ! local potential
2639 12 : nloc = sgppot%n_local
2640 48 : ALLOCATE (aloc(nloc), cloc(nloc))
2641 156 : aloc(1:nloc) = sgppot%a_local(1:nloc)
2642 156 : cloc(1:nloc) = sgppot%c_local(1:nloc)
2643 12 : CALL set_potential(qs_kind%sgp_potential, n_local=nloc, a_local=aloc, c_local=cloc)
2644 12 : IF (sgppot%has_nlcc) THEN
2645 0 : nlcc = sgppot%n_nlcc
2646 0 : ALLOCATE (anlcc(nlcc), cnlcc(nlcc))
2647 0 : anlcc(1:nlcc) = sgppot%a_nlcc(1:nlcc)
2648 0 : cnlcc(1:nlcc) = sgppot%c_nlcc(1:nlcc)
2649 0 : CALL set_potential(qs_kind%sgp_potential, has_nlcc=.TRUE., n_nlcc=nlcc, a_nlcc=anlcc, c_nlcc=cnlcc)
2650 : END IF
2651 12 : CALL set_potential(qs_kind%sgp_potential, z=z)
2652 12 : CALL get_qs_kind(qs_kind, elec_conf=elec_conf)
2653 12 : IF (.NOT. ASSOCIATED(elec_conf)) THEN
2654 12 : CALL set_qs_kind(qs_kind, elec_conf=upfpot%econf)
2655 : END IF
2656 12 : CALL get_qs_kind(qs_kind, elec_conf=elec_conf)
2657 12 : CALL set_potential(qs_kind%sgp_potential, elec_conf=elec_conf)
2658 12 : CALL atom_release_upf(upfpot)
2659 12 : CALL atom_sgp_release(sgppot)
2660 : CASE ("CNEO")
2661 8 : IF (zeff_correction /= 0.0_dp) THEN
2662 0 : CPABORT("CORE_CORRECTION is not compatible with CNEO")
2663 : END IF
2664 8 : CALL allocate_cneo_potential(qs_kind%cneo_potential)
2665 8 : CALL set_cneo_potential(qs_kind%cneo_potential, z=z)
2666 8 : mass = 0.0_dp
2667 : ! Input mass is the mass of the neutral atom, not the nucleus.
2668 : ! The mass of electrons will get subtracted later.
2669 : ! In principle, the most abundant pure isotope mass should be used.
2670 8 : IF (k_rep > 0) THEN
2671 : CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
2672 8 : keyword_name="MASS", n_rep_val=i)
2673 8 : IF (i > 0) CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
2674 2 : keyword_name="MASS", r_val=mass)
2675 : END IF
2676 : ! Remove electron mass from atomic mass to get nuclear mass
2677 8 : IF (mass - REAL(z, dp)*0.000548579909_dp > 0.0_dp) THEN
2678 2 : mass = mass - REAL(z, dp)*0.000548579909_dp
2679 2 : CALL set_cneo_potential(qs_kind%cneo_potential, mass=mass)
2680 : END IF
2681 : ! In case the mass is not set by user, get the default mass from z
2682 8 : CALL get_cneo_potential(qs_kind%cneo_potential, mass=mass)
2683 : ! Warn if the mass is from ptable
2684 8 : IF (ABS(mass + REAL(z, dp)*0.000548579909_dp - ptable(z)%amass) < 1.e-4_dp) THEN
2685 : CALL cp_warn(__LOCATION__, &
2686 : "Atomic mass of the atomic kind <"//TRIM(qs_kind%name)// &
2687 : "> is very close to its average mass. Is it a pure isotope? "// &
2688 : "Pure isotopes are preferable for CNEO. "// &
2689 0 : "(e.g., mass of 1H is 1.007825, not 1.00794)")
2690 : END IF
2691 8 : CALL get_qs_kind(qs_kind, elec_conf=elec_conf)
2692 8 : IF (.NOT. ASSOCIATED(elec_conf)) THEN
2693 8 : CALL get_cneo_potential(potential=qs_kind%cneo_potential, elec_conf=elec_conf)
2694 8 : CALL set_qs_kind(qs_kind, elec_conf=elec_conf)
2695 : ELSE
2696 0 : CALL set_cneo_potential(potential=qs_kind%cneo_potential, elec_conf=elec_conf)
2697 : END IF
2698 : CASE DEFAULT
2699 : CALL cp_abort(__LOCATION__, &
2700 : "An invalid potential type <"// &
2701 : TRIM(potential_name)//"> was specified "// &
2702 : "for the atomic kind <"// &
2703 10855 : TRIM(qs_kind%name))
2704 : END SELECT
2705 : ELSE
2706 : CALL cp_abort(__LOCATION__, &
2707 : "No potential type was defined for the "// &
2708 0 : "atomic kind <"//TRIM(qs_kind%name)//">")
2709 : END IF
2710 :
2711 10847 : CALL check_potential_basis_compatibility(qs_kind)
2712 :
2713 : ! Allocate and initialise the potential data set structure
2714 10847 : IF ((kgpot_name /= '') .OR. explicit_kgpot) THEN
2715 10847 : ipos = INDEX(kgpot_name, "-")
2716 10847 : IF (ipos > 1) THEN
2717 28 : kgpot_type = kgpot_name(:ipos - 1)
2718 : ELSE
2719 10819 : kgpot_type = kgpot_name
2720 : END IF
2721 10847 : CALL uppercase(kgpot_type)
2722 :
2723 10875 : SELECT CASE (TRIM(kgpot_type))
2724 : CASE ("TNADD")
2725 : ! determine the pseudopotential file to search
2726 28 : IF (kg_potential_fn_kind == "-") THEN
2727 28 : CALL section_vals_val_get(dft_section, "POTENTIAL_FILE_NAME", c_val=potential_file_name)
2728 : ELSE
2729 0 : potential_file_name = kg_potential_fn_kind
2730 : END IF
2731 28 : CALL allocate_potential(qs_kind%tnadd_potential)
2732 : CALL read_potential(qs_kind%element_symbol, kgpot_name, &
2733 : qs_kind%tnadd_potential, para_env, &
2734 28 : potential_file_name, kgpot_section, update_input)
2735 : CASE ("NONE")
2736 10819 : NULLIFY (qs_kind%tnadd_potential)
2737 : CASE DEFAULT
2738 : CALL cp_abort(__LOCATION__, &
2739 : "An invalid kg_potential type <"// &
2740 : TRIM(potential_name)//"> was specified "// &
2741 : "for the atomic kind <"// &
2742 10847 : TRIM(qs_kind%name))
2743 : END SELECT
2744 : END IF
2745 : END IF
2746 :
2747 : ! check that we have a nuclear orbital basis set if CNEO is requested
2748 10997 : nobasis_nuc = ASSOCIATED(qs_kind%cneo_potential)
2749 230937 : DO i = 1, SIZE(qs_kind%basis_sets)
2750 219940 : NULLIFY (tmp_basis_set)
2751 : CALL get_basis_from_container(qs_kind%basis_sets, basis_set=tmp_basis_set, &
2752 219940 : inumbas=i, basis_type=basis_type)
2753 230937 : IF (basis_type == "NUC") THEN
2754 8 : nobasis_nuc = .FALSE.
2755 8 : IF (.NOT. ASSOCIATED(qs_kind%cneo_potential)) THEN
2756 : CALL cp_warn(__LOCATION__, &
2757 : "POTENTIAL is not set to CNEO, NUC type basis set for KIND <"// &
2758 0 : TRIM(qs_kind%name)//"> will be ignored!")
2759 : END IF
2760 : END IF
2761 : END DO
2762 29688 : IF (nobasis_nuc) THEN
2763 : CALL cp_abort(__LOCATION__, &
2764 : "No NUC type basis set was defined for the "// &
2765 : "atomic kind <"//TRIM(qs_kind%name)//">, which is required by "// &
2766 0 : "POTENTIAL CNEO.")
2767 : END IF
2768 : END SELECT
2769 :
2770 16447 : CALL timestop(handle)
2771 :
2772 9835306 : END SUBROUTINE read_qs_kind
2773 :
2774 : ! **************************************************************************************************
2775 : !> \brief Ensure pseudo-potential and basis set were optimized for same number of valence electrons
2776 : !> \param qs_kind ...
2777 : !> \author Ole Schuett
2778 : ! **************************************************************************************************
2779 10847 : SUBROUTINE check_potential_basis_compatibility(qs_kind)
2780 : TYPE(qs_kind_type), INTENT(INOUT) :: qs_kind
2781 :
2782 : CHARACTER(LEN=default_string_length) :: name
2783 : INTEGER :: nbs, npp
2784 : TYPE(gth_potential_type), POINTER :: gth_potential
2785 : TYPE(gto_basis_set_type), POINTER :: basis_set
2786 :
2787 10847 : CALL get_qs_kind(qs_kind, name=name, gth_potential=gth_potential, basis_set=basis_set)
2788 :
2789 10847 : npp = -1; nbs = -1
2790 10847 : IF (ASSOCIATED(gth_potential)) THEN
2791 9461 : npp = parse_valence_electrons(gth_potential%aliases)
2792 : END IF
2793 10847 : IF (ASSOCIATED(basis_set)) THEN
2794 10847 : nbs = parse_valence_electrons(basis_set%aliases)
2795 : END IF
2796 :
2797 10847 : IF (npp >= 0 .AND. nbs >= 0 .AND. npp /= nbs) THEN
2798 : CALL cp_abort(__LOCATION__, "Basis-set and pseudo-potential of atomic kind '"//TRIM(name)//"'"// &
2799 0 : " were optimized for different valence electron numbers.")
2800 : END IF
2801 :
2802 10847 : END SUBROUTINE check_potential_basis_compatibility
2803 :
2804 : ! **************************************************************************************************
2805 : !> \brief Tries to parse valence eletron number using "-QXXX" notation, returns -1 if not found.
2806 : !> \param string ...
2807 : !> \return ...
2808 : !> \author Ole Schuett
2809 : ! **************************************************************************************************
2810 20308 : FUNCTION parse_valence_electrons(string) RESULT(n)
2811 : CHARACTER(*) :: string
2812 : INTEGER :: n
2813 :
2814 : INTEGER :: i, istat, j
2815 :
2816 20308 : i = INDEX(string, "-Q", .TRUE.)
2817 20308 : IF (i == 0) THEN
2818 6982 : n = -1
2819 : ELSE
2820 13326 : j = SCAN(string(i + 2:), "- ")
2821 13326 : READ (string(i + 2:i + j), '(I3)', iostat=istat) n
2822 13326 : IF (istat /= 0) n = -1
2823 : END IF
2824 :
2825 20308 : END FUNCTION parse_valence_electrons
2826 :
2827 : ! **************************************************************************************************
2828 : !> \brief Read an atomic kind set data set from the input file.
2829 : !> \param qs_kind_set ...
2830 : !> \param atomic_kind_set ...
2831 : !> \param kind_section ...
2832 : !> \param para_env ...
2833 : !> \param force_env_section ...
2834 : !> \param silent ...
2835 : ! **************************************************************************************************
2836 8712 : SUBROUTINE create_qs_kind_set(qs_kind_set, atomic_kind_set, kind_section, para_env, &
2837 : force_env_section, silent)
2838 :
2839 : TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
2840 : TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
2841 : TYPE(section_vals_type), POINTER :: kind_section
2842 : TYPE(mp_para_env_type), POINTER :: para_env
2843 : TYPE(section_vals_type), POINTER :: force_env_section
2844 : LOGICAL, INTENT(IN) :: silent
2845 :
2846 : CHARACTER(len=*), PARAMETER :: routineN = 'create_qs_kind_set'
2847 :
2848 : INTEGER :: handle, ikind, method, nkind, qs_method
2849 : LOGICAL :: no_fail
2850 :
2851 8712 : CALL timeset(routineN, handle)
2852 :
2853 8712 : IF (ASSOCIATED(qs_kind_set)) CPABORT("create_qs_kind_set: qs_kind_set already associated")
2854 8712 : IF (.NOT. ASSOCIATED(atomic_kind_set)) CPABORT("create_qs_kind_set: atomic_kind_set not associated")
2855 :
2856 8712 : no_fail = .FALSE.
2857 :
2858 : ! Between all methods only SE and DFTB/xTB may not need a KIND section.
2859 8712 : CALL section_vals_val_get(force_env_section, "METHOD", i_val=method)
2860 8712 : IF (method == do_qs) THEN
2861 8688 : CALL section_vals_val_get(force_env_section, "DFT%QS%METHOD", i_val=qs_method)
2862 1000 : SELECT CASE (qs_method)
2863 : CASE (do_method_mndo, do_method_am1, do_method_pm3, do_method_pm6fm, do_method_pm6, &
2864 : do_method_pdg, do_method_rm1, do_method_mndod, do_method_pnnl)
2865 1000 : no_fail = .TRUE.
2866 : CASE (do_method_dftb)
2867 292 : no_fail = .TRUE.
2868 : CASE (do_method_xtb)
2869 8688 : no_fail = .TRUE.
2870 : END SELECT
2871 24 : ELSE IF (method == do_sirius) THEN
2872 20 : qs_method = do_method_pw
2873 : ELSE
2874 4 : qs_method = method
2875 : END IF
2876 :
2877 8712 : nkind = SIZE(atomic_kind_set)
2878 216823 : ALLOCATE (qs_kind_set(nkind))
2879 :
2880 25159 : DO ikind = 1, nkind
2881 16447 : qs_kind_set(ikind)%name = atomic_kind_set(ikind)%name
2882 16447 : qs_kind_set(ikind)%element_symbol = atomic_kind_set(ikind)%element_symbol
2883 16447 : qs_kind_set(ikind)%natom = atomic_kind_set(ikind)%natom
2884 : CALL read_qs_kind(qs_kind_set(ikind), kind_section, para_env, force_env_section, &
2885 25159 : no_fail, qs_method, silent)
2886 : END DO
2887 :
2888 8712 : CALL timestop(handle)
2889 :
2890 17424 : END SUBROUTINE create_qs_kind_set
2891 :
2892 : ! **************************************************************************************************
2893 : !> \brief This routines should perform only checks. no settings are allowed at
2894 : !> this level anymore..
2895 : !> \param qs_kind ...
2896 : !> \param dft_control ...
2897 : !> \param subsys_section ...
2898 : ! **************************************************************************************************
2899 16365 : SUBROUTINE check_qs_kind(qs_kind, dft_control, subsys_section)
2900 :
2901 : TYPE(qs_kind_type), POINTER :: qs_kind
2902 : TYPE(dft_control_type), INTENT(IN) :: dft_control
2903 : TYPE(section_vals_type), POINTER :: subsys_section
2904 :
2905 : INTEGER :: gfn_type
2906 : LOGICAL :: defined
2907 : TYPE(qs_dftb_atom_type), POINTER :: dftb_parameter
2908 : TYPE(semi_empirical_type), POINTER :: se_parameter
2909 : TYPE(xtb_atom_type), POINTER :: xtb_parameter
2910 :
2911 16365 : IF (dft_control%qs_control%semi_empirical) THEN
2912 2244 : CALL get_qs_kind(qs_kind, se_parameter=se_parameter)
2913 2244 : CPASSERT(ASSOCIATED(se_parameter))
2914 2244 : CALL get_se_param(se_parameter, defined=defined)
2915 2244 : CPASSERT(defined)
2916 2244 : CALL write_se_param(se_parameter, subsys_section)
2917 14121 : ELSE IF (dft_control%qs_control%dftb) THEN
2918 618 : CALL get_qs_kind(qs_kind, dftb_parameter=dftb_parameter)
2919 618 : CPASSERT(ASSOCIATED(dftb_parameter))
2920 618 : CALL get_dftb_atom_param(dftb_parameter, defined=defined)
2921 618 : CPASSERT(defined)
2922 618 : CALL write_dftb_atom_param(dftb_parameter, subsys_section)
2923 13503 : ELSE IF (dft_control%qs_control%xtb) THEN
2924 2562 : IF (.NOT. (dft_control%qs_control%xtb_control%do_tblite)) THEN
2925 2256 : CALL get_qs_kind(qs_kind, xtb_parameter=xtb_parameter)
2926 2256 : CPASSERT(ASSOCIATED(xtb_parameter))
2927 2256 : gfn_type = dft_control%qs_control%xtb_control%gfn_type
2928 2256 : CALL write_xtb_atom_param(xtb_parameter, gfn_type, subsys_section)
2929 : END IF
2930 : END IF
2931 :
2932 16365 : END SUBROUTINE check_qs_kind
2933 :
2934 : ! **************************************************************************************************
2935 : !> \brief ...
2936 : !> \param qs_kind_set ...
2937 : !> \param dft_control ...
2938 : !> \param subsys_section ...
2939 : ! **************************************************************************************************
2940 8664 : SUBROUTINE check_qs_kind_set(qs_kind_set, dft_control, subsys_section)
2941 :
2942 : TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
2943 : TYPE(dft_control_type), INTENT(IN) :: dft_control
2944 : TYPE(section_vals_type), POINTER :: subsys_section
2945 :
2946 : CHARACTER(len=*), PARAMETER :: routineN = 'check_qs_kind_set'
2947 :
2948 : INTEGER :: handle, ikind, nkind
2949 : TYPE(qs_kind_type), POINTER :: qs_kind
2950 :
2951 8664 : CALL timeset(routineN, handle)
2952 8664 : IF (ASSOCIATED(qs_kind_set)) THEN
2953 8664 : nkind = SIZE(qs_kind_set)
2954 25029 : DO ikind = 1, nkind
2955 16365 : qs_kind => qs_kind_set(ikind)
2956 25029 : CALL check_qs_kind(qs_kind, dft_control, subsys_section)
2957 : END DO
2958 8664 : IF (dft_control%qs_control%xtb) THEN
2959 : CALL write_xtb_kab_param(qs_kind_set, subsys_section, &
2960 1192 : dft_control%qs_control%xtb_control)
2961 : END IF
2962 : ELSE
2963 0 : CPABORT("The pointer qs_kind_set is not associated")
2964 : END IF
2965 8664 : CALL timestop(handle)
2966 8664 : END SUBROUTINE check_qs_kind_set
2967 :
2968 : ! **************************************************************************************************
2969 : !> \brief ...
2970 : !> \param qs_kind_set ...
2971 : !> \param subsys_section ...
2972 : !> \param xtb_control ...
2973 : ! **************************************************************************************************
2974 1192 : SUBROUTINE write_xtb_kab_param(qs_kind_set, subsys_section, xtb_control)
2975 :
2976 : TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
2977 : TYPE(section_vals_type), POINTER :: subsys_section
2978 : TYPE(xtb_control_type), POINTER :: xtb_control
2979 :
2980 : CHARACTER(LEN=default_string_length) :: aname, bname
2981 : INTEGER :: ikind, io_unit, jkind, nkind, za, zb
2982 : TYPE(cp_logger_type), POINTER :: logger
2983 : TYPE(qs_kind_type), POINTER :: qs_kinda, qs_kindb
2984 : TYPE(xtb_atom_type), POINTER :: xtb_parameter_a, xtb_parameter_b
2985 :
2986 1192 : NULLIFY (logger)
2987 1192 : logger => cp_get_default_logger()
2988 1192 : IF (BTEST(cp_print_key_should_output(logger%iter_info, subsys_section, &
2989 : "PRINT%KINDS/POTENTIAL"), cp_p_file)) THEN
2990 :
2991 0 : io_unit = cp_print_key_unit_nr(logger, subsys_section, "PRINT%KINDS", extension=".Log")
2992 0 : IF (io_unit > 0) THEN
2993 :
2994 0 : WRITE (io_unit, "(/,T2,A)") "xTB| Kab parameters"
2995 0 : nkind = SIZE(qs_kind_set)
2996 0 : DO ikind = 1, nkind
2997 0 : qs_kinda => qs_kind_set(ikind)
2998 0 : CALL get_qs_kind(qs_kinda, xtb_parameter=xtb_parameter_a)
2999 0 : CALL get_xtb_atom_param(xtb_parameter_a, aname=aname, z=za)
3000 0 : DO jkind = ikind, nkind
3001 0 : qs_kindb => qs_kind_set(jkind)
3002 0 : CALL get_qs_kind(qs_kindb, xtb_parameter=xtb_parameter_b)
3003 0 : CALL get_xtb_atom_param(xtb_parameter_b, aname=bname, z=zb)
3004 : WRITE (io_unit, "(A,T10,A15,T25,A15,T71,F10.3)") &
3005 0 : " Kab:", TRIM(aname), TRIM(bname), xtb_set_kab(za, zb, xtb_control)
3006 : END DO
3007 : END DO
3008 0 : WRITE (io_unit, *)
3009 :
3010 : END IF
3011 :
3012 0 : CALL cp_print_key_finished_output(io_unit, logger, subsys_section, "PRINT%KINDS")
3013 : END IF
3014 :
3015 1192 : END SUBROUTINE write_xtb_kab_param
3016 :
3017 : ! **************************************************************************************************
3018 : !> \brief Set the components of an atomic kind data set.
3019 : !> \param qs_kind ...
3020 : !> \param paw_atom ...
3021 : !> \param ghost ...
3022 : !> \param floating ...
3023 : !> \param hard_radius ...
3024 : !> \param hard0_radius ...
3025 : !> \param covalent_radius ...
3026 : !> \param vdw_radius ...
3027 : !> \param lmax_rho0 ...
3028 : !> \param zeff ...
3029 : !> \param no_optimize ...
3030 : !> \param dispersion ...
3031 : !> \param u_minus_j ...
3032 : !> \param reltmat ...
3033 : !> \param dftb_parameter ...
3034 : !> \param xtb_parameter ...
3035 : !> \param elec_conf ...
3036 : !> \param pao_basis_size ...
3037 : ! **************************************************************************************************
3038 30199 : SUBROUTINE set_qs_kind(qs_kind, paw_atom, ghost, floating, hard_radius, hard0_radius, &
3039 : covalent_radius, vdw_radius, lmax_rho0, zeff, &
3040 : no_optimize, dispersion, u_minus_j, reltmat, &
3041 : dftb_parameter, xtb_parameter, &
3042 30199 : elec_conf, pao_basis_size)
3043 :
3044 : TYPE(qs_kind_type), INTENT(INOUT) :: qs_kind
3045 : LOGICAL, INTENT(IN), OPTIONAL :: paw_atom, ghost, floating
3046 : REAL(KIND=dp), INTENT(IN), OPTIONAL :: hard_radius, hard0_radius, &
3047 : covalent_radius, vdw_radius
3048 : INTEGER, INTENT(IN), OPTIONAL :: lmax_rho0
3049 : REAL(KIND=dp), INTENT(IN), OPTIONAL :: zeff
3050 : LOGICAL, INTENT(IN), OPTIONAL :: no_optimize
3051 : TYPE(qs_atom_dispersion_type), OPTIONAL, POINTER :: dispersion
3052 : REAL(KIND=dp), INTENT(IN), OPTIONAL :: u_minus_j
3053 : REAL(KIND=dp), DIMENSION(:, :), OPTIONAL, POINTER :: reltmat
3054 : TYPE(qs_dftb_atom_type), OPTIONAL, POINTER :: dftb_parameter
3055 : TYPE(xtb_atom_type), OPTIONAL, POINTER :: xtb_parameter
3056 : INTEGER, DIMENSION(:), INTENT(IN), OPTIONAL :: elec_conf
3057 : INTEGER, INTENT(IN), OPTIONAL :: pao_basis_size
3058 :
3059 30199 : IF (PRESENT(dftb_parameter)) qs_kind%dftb_parameter => dftb_parameter
3060 30199 : IF (PRESENT(xtb_parameter)) qs_kind%xtb_parameter => xtb_parameter
3061 30199 : IF (PRESENT(elec_conf)) THEN
3062 16297 : IF (ASSOCIATED(qs_kind%elec_conf)) THEN
3063 0 : DEALLOCATE (qs_kind%elec_conf)
3064 : END IF
3065 48891 : ALLOCATE (qs_kind%elec_conf(0:SIZE(elec_conf) - 1))
3066 59619 : qs_kind%elec_conf(:) = elec_conf(:)
3067 : END IF
3068 30199 : IF (PRESENT(paw_atom)) qs_kind%paw_atom = paw_atom
3069 30199 : IF (PRESENT(hard_radius)) qs_kind%hard_radius = hard_radius
3070 30199 : IF (PRESENT(hard0_radius)) qs_kind%hard0_radius = hard0_radius
3071 30199 : IF (PRESENT(covalent_radius)) qs_kind%covalent_radius = covalent_radius
3072 30199 : IF (PRESENT(vdw_radius)) qs_kind%vdw_radius = vdw_radius
3073 30199 : IF (PRESENT(lmax_rho0)) qs_kind%lmax_rho0 = lmax_rho0
3074 30199 : IF (PRESENT(zeff)) THEN
3075 0 : IF (ASSOCIATED(qs_kind%all_potential)) THEN
3076 0 : CALL set_potential(potential=qs_kind%all_potential, zeff=zeff)
3077 0 : ELSE IF (ASSOCIATED(qs_kind%gth_potential)) THEN
3078 0 : CALL set_potential(potential=qs_kind%gth_potential, zeff=zeff)
3079 0 : ELSE IF (ASSOCIATED(qs_kind%sgp_potential)) THEN
3080 0 : CALL set_potential(potential=qs_kind%sgp_potential, zeff=zeff)
3081 0 : ELSE IF (ASSOCIATED(qs_kind%cneo_potential)) THEN
3082 0 : CPABORT("CNEO potential ZEFF should not be manually set")
3083 : END IF
3084 : END IF
3085 30199 : IF (PRESENT(ghost)) qs_kind%ghost = ghost
3086 :
3087 30199 : IF (PRESENT(floating)) qs_kind%floating = floating
3088 :
3089 30199 : IF (PRESENT(no_optimize)) qs_kind%no_optimize = no_optimize
3090 :
3091 30199 : IF (PRESENT(dispersion)) qs_kind%dispersion => dispersion
3092 :
3093 30199 : IF (PRESENT(u_minus_j)) THEN
3094 476 : IF (ASSOCIATED(qs_kind%dft_plus_u)) THEN
3095 476 : qs_kind%dft_plus_u%u_minus_j = u_minus_j
3096 : END IF
3097 : END IF
3098 :
3099 30199 : IF (PRESENT(reltmat)) qs_kind%reltmat => reltmat
3100 :
3101 30199 : IF (PRESENT(pao_basis_size)) qs_kind%pao_basis_size = pao_basis_size
3102 :
3103 30199 : END SUBROUTINE set_qs_kind
3104 :
3105 : ! **************************************************************************************************
3106 : !> \brief Write an atomic kind data set to the output unit.
3107 : !> \param qs_kind ...
3108 : !> \param kind_number ...
3109 : !> \param output_unit ...
3110 : !> \par History
3111 : !> Creation (09.02.2002,MK)
3112 : ! **************************************************************************************************
3113 4010 : SUBROUTINE write_qs_kind(qs_kind, kind_number, output_unit)
3114 :
3115 : TYPE(qs_kind_type), POINTER :: qs_kind
3116 : INTEGER, INTENT(in) :: kind_number, output_unit
3117 :
3118 : CHARACTER(LEN=3) :: yon
3119 : CHARACTER(LEN=default_string_length) :: basis_type, bstring
3120 : INTEGER :: ibas
3121 : LOGICAL :: do_print
3122 : TYPE(gto_basis_set_type), POINTER :: tmp_basis
3123 :
3124 4010 : IF (output_unit > 0) THEN
3125 :
3126 4010 : IF (ASSOCIATED(qs_kind)) THEN
3127 : WRITE (UNIT=output_unit, FMT="(/,T2,I2,A,T57,A,T75,I6)") &
3128 4010 : kind_number, ". Atomic kind: "//TRIM(qs_kind%name), &
3129 8020 : "Number of atoms: ", qs_kind%natom
3130 :
3131 84210 : DO ibas = 1, SIZE(qs_kind%basis_sets, 1)
3132 80200 : NULLIFY (tmp_basis)
3133 : CALL get_basis_from_container(qs_kind%basis_sets, basis_set=tmp_basis, &
3134 80200 : inumbas=ibas, basis_type=basis_type)
3135 80200 : do_print = .TRUE.
3136 74994 : SELECT CASE (basis_type)
3137 : CASE DEFAULT
3138 74994 : bstring = "Basis Set"
3139 3907 : do_print = .FALSE.
3140 : CASE ("ORB")
3141 3907 : bstring = "Orbital Basis Set"
3142 : CASE ("ORB_SOFT")
3143 583 : bstring = "GAPW Soft Basis Set"
3144 0 : do_print = .FALSE.
3145 : CASE ("AUX")
3146 0 : bstring = "Auxiliary Basis Set"
3147 : CASE ("MIN")
3148 0 : bstring = "Minimal Basis Set"
3149 : CASE ("RI_AUX")
3150 352 : bstring = "RI Auxiliary Basis Set"
3151 : CASE ("AUX_FIT")
3152 244 : bstring = "Auxiliary Fit Basis Set"
3153 : CASE ("LRI_AUX")
3154 23 : bstring = "LRI Basis Set"
3155 : CASE ("P_LRI_AUX")
3156 4 : bstring = "LRI Basis Set for TDDFPT"
3157 : CASE ("RI_XAS")
3158 0 : bstring = "RI XAS Basis Set"
3159 : CASE ("RI_HFX")
3160 85 : bstring = "RI HFX Basis Set"
3161 : CASE ("NUC")
3162 4 : bstring = "Nuclear Basis Set"
3163 4 : do_print = .FALSE.
3164 : CASE ("NUC_SOFT")
3165 4 : bstring = "Nuclear Soft Basis Set"
3166 80200 : do_print = .FALSE.
3167 : END SELECT
3168 :
3169 4010 : IF (do_print) THEN
3170 4615 : CALL write_orb_basis_set(tmp_basis, output_unit, bstring)
3171 : END IF
3172 :
3173 : END DO
3174 :
3175 4010 : IF (qs_kind%ghost) THEN
3176 : WRITE (UNIT=output_unit, FMT="(/,T6,A)") &
3177 11 : "The atoms of this atomic kind are GHOST atoms!"
3178 : END IF
3179 4010 : IF (qs_kind%monovalent) THEN
3180 : WRITE (UNIT=output_unit, FMT="(/,T6,A)") &
3181 1 : "The atoms of this atomic kind are MONOVALENT!"
3182 : END IF
3183 4010 : IF (qs_kind%floating) THEN
3184 : WRITE (UNIT=output_unit, FMT="(/,T6,A)") &
3185 0 : "The atoms of this atomic kind are FLOATING BASIS FUNCTIONS."
3186 : END IF
3187 4010 : IF (qs_kind%covalent_radius > 0.0_dp) THEN
3188 : WRITE (UNIT=output_unit, FMT="(/,T8,A,T71,F10.3)") &
3189 2713 : "Atomic covalent radius [Angstrom]:", &
3190 5426 : qs_kind%covalent_radius*angstrom
3191 : END IF
3192 4010 : IF (qs_kind%vdw_radius > 0.0_dp) THEN
3193 : WRITE (UNIT=output_unit, FMT="(/,T8,A,T71,F10.3)") &
3194 2713 : "Atomic van der Waals radius [Angstrom]:", &
3195 5426 : qs_kind%vdw_radius*angstrom
3196 : END IF
3197 4010 : IF (qs_kind%paw_atom) THEN
3198 : WRITE (UNIT=output_unit, FMT="(/,T6,A)") &
3199 492 : "The atoms of this atomic kind are PAW atoms (GAPW):"
3200 : WRITE (UNIT=output_unit, FMT="(T8,A,T71,F10.3)") &
3201 492 : "Hard Gaussian function radius:", qs_kind%hard_radius, &
3202 492 : "Rho0 radius:", qs_kind%hard0_radius, &
3203 492 : "Maximum GTO radius used for PAW projector construction:", &
3204 984 : qs_kind%max_rad_local
3205 492 : NULLIFY (tmp_basis)
3206 : CALL get_basis_from_container(qs_kind%basis_sets, basis_set=tmp_basis, &
3207 492 : basis_type="ORB_SOFT")
3208 492 : CALL write_orb_basis_set(tmp_basis, output_unit, "GAPW Soft Basis Set")
3209 : END IF
3210 : ! Potentials
3211 4010 : IF (ASSOCIATED(qs_kind%all_potential)) CALL write_potential(qs_kind%all_potential, output_unit)
3212 4010 : IF (ASSOCIATED(qs_kind%gth_potential)) CALL write_potential(qs_kind%gth_potential, output_unit)
3213 4010 : IF (ASSOCIATED(qs_kind%sgp_potential)) CALL write_potential(qs_kind%sgp_potential, output_unit)
3214 4010 : IF (ASSOCIATED(qs_kind%tnadd_potential)) CALL write_potential(qs_kind%tnadd_potential, output_unit)
3215 4010 : IF (ASSOCIATED(qs_kind%dft_plus_u)) THEN
3216 : WRITE (UNIT=output_unit, FMT="(/,T6,A,/,T8,A,T76,I5,/,T8,A,T73,F8.3)") &
3217 16 : "A DFT+U correction is applied to atoms of this atomic kind:", &
3218 16 : "Angular quantum momentum number L:", qs_kind%dft_plus_u%l, &
3219 32 : "U(eff) = (U - J) value in [eV]:", qs_kind%dft_plus_u%u_minus_j_target*evolt
3220 16 : IF (qs_kind%dft_plus_u%u_ramping > 0.0_dp) THEN
3221 4 : IF (qs_kind%dft_plus_u%init_u_ramping_each_scf) THEN
3222 2 : yon = "YES"
3223 : ELSE
3224 2 : yon = " NO"
3225 : END IF
3226 : WRITE (UNIT=output_unit, FMT="(T8,A,T73,F8.3,/,T8,A,T73,ES8.1,/,T8,A,T78,A3)") &
3227 4 : "Increment for U ramping in [eV]:", qs_kind%dft_plus_u%u_ramping*evolt, &
3228 4 : "SCF threshold value for U ramping:", qs_kind%dft_plus_u%eps_u_ramping, &
3229 8 : "Set U ramping value to zero before each wavefunction optimisation:", yon
3230 : END IF
3231 16 : IF (ASSOCIATED(qs_kind%dft_plus_u%orbitals)) THEN
3232 : WRITE (UNIT=output_unit, FMT="(T8,A)") &
3233 2 : "An initial orbital occupation is requested:"
3234 2 : IF (ASSOCIATED(qs_kind%dft_plus_u%nelec)) THEN
3235 4 : IF (ANY(qs_kind%dft_plus_u%nelec(:) >= 0.5_dp)) THEN
3236 0 : IF (SIZE(qs_kind%dft_plus_u%nelec) > 1) THEN
3237 : WRITE (UNIT=output_unit, FMT="(T9,A,T75,F6.2)") &
3238 0 : "Number of alpha electrons:", &
3239 0 : qs_kind%dft_plus_u%nelec(1), &
3240 0 : "Number of beta electrons:", &
3241 0 : qs_kind%dft_plus_u%nelec(2)
3242 : ELSE
3243 : WRITE (UNIT=output_unit, FMT="(T9,A,T75,F6.2)") &
3244 0 : "Number of electrons:", &
3245 0 : qs_kind%dft_plus_u%nelec(1)
3246 : END IF
3247 : END IF
3248 : END IF
3249 : WRITE (UNIT=output_unit, FMT="(T9,A,(T78,I3))") &
3250 2 : "Preferred (initial) orbital occupation order (orbital M values):", &
3251 10 : qs_kind%dft_plus_u%orbitals(:)
3252 : WRITE (UNIT=output_unit, FMT="(T9,A,T71,ES10.3,/,T9,A,T76,I5)") &
3253 2 : "Threshold value for the SCF convergence criterion:", &
3254 2 : qs_kind%dft_plus_u%eps_scf, &
3255 2 : "Number of initial SCF iterations:", &
3256 4 : qs_kind%dft_plus_u%max_scf
3257 2 : IF (qs_kind%dft_plus_u%smear) THEN
3258 : WRITE (UNIT=output_unit, FMT="(T9,A)") &
3259 2 : "A smearing of the orbital occupations will be performed"
3260 : END IF
3261 : END IF
3262 : END IF
3263 4010 : IF (ASSOCIATED(qs_kind%cneo_potential)) THEN
3264 : WRITE (UNIT=output_unit, FMT="(/,T6,A)") &
3265 4 : "The nuclei of this atomic kind are quantum mechanical (CNEO)"
3266 4 : CALL write_cneo_potential(qs_kind%cneo_potential, output_unit)
3267 4 : NULLIFY (tmp_basis)
3268 : CALL get_basis_from_container(qs_kind%basis_sets, basis_set=tmp_basis, &
3269 4 : basis_type="NUC")
3270 4 : CALL write_orb_basis_set(tmp_basis, output_unit, "Nuclear Basis Set")
3271 4 : NULLIFY (tmp_basis)
3272 : CALL get_basis_from_container(qs_kind%basis_sets, basis_set=tmp_basis, &
3273 4 : basis_type="NUC_SOFT")
3274 4 : CALL write_orb_basis_set(tmp_basis, output_unit, "Nuclear Soft Basis Set")
3275 : END IF
3276 : ELSE
3277 0 : CPABORT("write_qs_kind requires qs_kind to be associated")
3278 : END IF
3279 :
3280 : END IF
3281 :
3282 4010 : END SUBROUTINE write_qs_kind
3283 :
3284 : ! **************************************************************************************************
3285 : !> \brief Write an atomic kind set data set to the output unit.
3286 : !> \param qs_kind_set ...
3287 : !> \param subsys_section ...
3288 : !> \par History
3289 : !> Creation (09.02.2002,MK)
3290 : ! **************************************************************************************************
3291 8678 : SUBROUTINE write_qs_kind_set(qs_kind_set, subsys_section)
3292 : TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
3293 : TYPE(section_vals_type), POINTER :: subsys_section
3294 :
3295 : CHARACTER(len=*), PARAMETER :: routineN = 'write_qs_kind_set'
3296 :
3297 : INTEGER :: handle, ikind, nkind, output_unit
3298 : TYPE(cp_logger_type), POINTER :: logger
3299 : TYPE(qs_kind_type), POINTER :: qs_kind
3300 :
3301 8678 : CALL timeset(routineN, handle)
3302 :
3303 8678 : NULLIFY (logger)
3304 8678 : logger => cp_get_default_logger()
3305 : output_unit = cp_print_key_unit_nr(logger, subsys_section, &
3306 8678 : "PRINT%KINDS", extension=".Log")
3307 8678 : IF (output_unit > 0) THEN
3308 2155 : IF (ASSOCIATED(qs_kind_set)) THEN
3309 2155 : WRITE (UNIT=output_unit, FMT="(/,/,T2,A)") "ATOMIC KIND INFORMATION"
3310 2155 : nkind = SIZE(qs_kind_set)
3311 6165 : DO ikind = 1, nkind
3312 4010 : qs_kind => qs_kind_set(ikind)
3313 6165 : CALL write_qs_kind(qs_kind, ikind, output_unit)
3314 : END DO
3315 : ELSE
3316 0 : CPABORT("write_qs_kind_set requires qs_kind_set to be associated")
3317 : END IF
3318 : END IF
3319 :
3320 : CALL cp_print_key_finished_output(output_unit, logger, subsys_section, &
3321 8678 : "PRINT%KINDS")
3322 :
3323 8678 : CALL timestop(handle)
3324 :
3325 8678 : END SUBROUTINE write_qs_kind_set
3326 :
3327 : ! **************************************************************************************************
3328 : !> \brief Write all the GTO basis sets of an atomic kind set to the output
3329 : !> unit (for the printing of the unnormalized basis sets as read from
3330 : !> database).
3331 : !> \param qs_kind_set ...
3332 : !> \param subsys_section ...
3333 : !> \par History
3334 : !> Creation (17.01.2002,MK)
3335 : ! **************************************************************************************************
3336 8658 : SUBROUTINE write_gto_basis_sets(qs_kind_set, subsys_section)
3337 :
3338 : TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
3339 : TYPE(section_vals_type), POINTER :: subsys_section
3340 :
3341 : CHARACTER(LEN=*), PARAMETER :: routineN = 'write_gto_basis_sets'
3342 :
3343 : CHARACTER(LEN=default_string_length) :: basis_type, bstring
3344 : INTEGER :: handle, ibas, ikind, nkind, output_unit
3345 : TYPE(cp_logger_type), POINTER :: logger
3346 : TYPE(gto_basis_set_type), POINTER :: tmp_basis
3347 : TYPE(qs_kind_type), POINTER :: qs_kind
3348 :
3349 8658 : CALL timeset(routineN, handle)
3350 :
3351 8658 : NULLIFY (logger)
3352 8658 : logger => cp_get_default_logger()
3353 : output_unit = cp_print_key_unit_nr(logger, subsys_section, &
3354 : "PRINT%KINDS/BASIS_SET", &
3355 8658 : extension=".Log")
3356 8658 : IF (output_unit > 0) THEN
3357 63 : IF (ASSOCIATED(qs_kind_set)) THEN
3358 : WRITE (UNIT=output_unit, FMT="(/,/,T2,A)") &
3359 63 : "BASIS SET INFORMATION (Unnormalised Gaussian-type functions)"
3360 63 : nkind = SIZE(qs_kind_set)
3361 181 : DO ikind = 1, nkind
3362 118 : qs_kind => qs_kind_set(ikind)
3363 : WRITE (UNIT=output_unit, FMT="(/,T2,I2,A)") &
3364 118 : ikind, ". Atomic kind: "//TRIM(qs_kind%name)
3365 :
3366 2541 : DO ibas = 1, SIZE(qs_kind%basis_sets, 1)
3367 2360 : NULLIFY (tmp_basis)
3368 : CALL get_basis_from_container(qs_kind%basis_sets, basis_set=tmp_basis, &
3369 2360 : inumbas=ibas, basis_type=basis_type)
3370 2360 : IF (basis_type == "") CYCLE
3371 11 : SELECT CASE (basis_type)
3372 : CASE DEFAULT
3373 11 : bstring = "Basis Set"
3374 : CASE ("ORB")
3375 118 : bstring = "Orbital Basis Set"
3376 : CASE ("ORB_SOFT")
3377 11 : bstring = "GAPW Soft Basis Set"
3378 : CASE ("AUX")
3379 0 : bstring = "Auxiliary Basis Set"
3380 : CASE ("MIN")
3381 0 : bstring = "Minimal Basis Set"
3382 : CASE ("RI_AUX")
3383 0 : bstring = "RI Auxiliary Basis Set"
3384 : CASE ("AUX_FIT")
3385 0 : bstring = "Auxiliary Fit Basis Set"
3386 : CASE ("LRI_AUX")
3387 2 : bstring = "LRI Basis Set"
3388 : CASE ("P_LRI_AUX")
3389 0 : bstring = "LRI Basis Set for TDDFPT"
3390 : CASE ("RI_HFX")
3391 0 : bstring = "RI HFX Basis Set"
3392 : CASE ("NUC")
3393 0 : bstring = "Nuclear Basis Set"
3394 0 : IF (.NOT. ASSOCIATED(qs_kind%cneo_potential)) NULLIFY (tmp_basis)
3395 : CASE ("NUC_SOFT")
3396 0 : bstring = "Nuclear Soft Basis Set"
3397 142 : IF (.NOT. ASSOCIATED(qs_kind%cneo_potential)) NULLIFY (tmp_basis)
3398 : END SELECT
3399 :
3400 260 : IF (ASSOCIATED(tmp_basis)) CALL write_gto_basis_set(tmp_basis, output_unit, bstring)
3401 :
3402 : END DO
3403 :
3404 : END DO
3405 : ELSE
3406 0 : CPABORT("write_gto_basis_sets require qs_kind_set to be associated")
3407 : END IF
3408 : END IF
3409 :
3410 : CALL cp_print_key_finished_output(output_unit, logger, subsys_section, &
3411 8658 : "PRINT%KINDS/BASIS_SET")
3412 :
3413 8658 : CALL timestop(handle)
3414 :
3415 8658 : END SUBROUTINE write_gto_basis_sets
3416 :
3417 : ! **************************************************************************************************
3418 : !> \brief ...
3419 : !> \param atomic_kind ...
3420 : !> \param qs_kind ...
3421 : !> \param ncalc ...
3422 : !> \param ncore ...
3423 : !> \param nelem ...
3424 : !> \param edelta ...
3425 : ! **************************************************************************************************
3426 105542 : SUBROUTINE init_atom_electronic_state(atomic_kind, qs_kind, ncalc, ncore, nelem, edelta)
3427 :
3428 : TYPE(atomic_kind_type), INTENT(IN) :: atomic_kind
3429 : TYPE(qs_kind_type), INTENT(IN) :: qs_kind
3430 : INTEGER, DIMENSION(0:lmat, 10), INTENT(OUT) :: ncalc, ncore, nelem
3431 : REAL(KIND=dp), DIMENSION(0:lmat, 10, 2), &
3432 : INTENT(OUT) :: edelta
3433 :
3434 : INTEGER :: i, ii, is, l, ll, ne, nn, z
3435 52771 : INTEGER, DIMENSION(:), POINTER :: econf
3436 52771 : INTEGER, DIMENSION(:, :), POINTER :: addel, laddel, naddel
3437 : LOGICAL :: bs_occupation, monovalent
3438 : REAL(KIND=dp) :: dmag, magnetization
3439 : TYPE(gth_potential_type), POINTER :: gth_potential
3440 : TYPE(sgp_potential_type), POINTER :: sgp_potential
3441 :
3442 52771 : CALL get_atomic_kind(atomic_kind, z=z)
3443 52771 : NULLIFY (gth_potential)
3444 : CALL get_qs_kind(qs_kind, &
3445 : gth_potential=gth_potential, &
3446 : sgp_potential=sgp_potential, &
3447 : magnetization=magnetization, &
3448 : bs_occupation=bs_occupation, &
3449 : monovalent=monovalent, &
3450 52771 : addel=addel, laddel=laddel, naddel=naddel)
3451 :
3452 : ! electronic state
3453 52771 : nelem = 0
3454 52771 : ncore = 0
3455 52771 : ncalc = 0
3456 52771 : edelta = 0.0_dp
3457 52771 : IF (monovalent) THEN
3458 4 : ncalc(0, 1) = 1
3459 4 : nelem(0, 1) = 1
3460 52767 : ELSE IF (ASSOCIATED(gth_potential)) THEN
3461 29209 : CALL get_potential(gth_potential, elec_conf=econf)
3462 29209 : CALL set_pseudo_state(econf, z, ncalc, ncore, nelem)
3463 23558 : ELSE IF (ASSOCIATED(sgp_potential)) THEN
3464 190 : CALL get_potential(sgp_potential, elec_conf=econf)
3465 190 : CALL set_pseudo_state(econf, z, ncalc, ncore, nelem)
3466 : ELSE
3467 116840 : DO l = 0, MIN(lmat, UBOUND(ptable(z)%e_conv, 1))
3468 93472 : ll = 2*(2*l + 1)
3469 93472 : nn = ptable(z)%e_conv(l)
3470 93472 : ii = 0
3471 23368 : DO
3472 132272 : ii = ii + 1
3473 132272 : IF (nn <= ll) THEN
3474 93472 : nelem(l, ii) = nn
3475 : EXIT
3476 : ELSE
3477 38800 : nelem(l, ii) = ll
3478 38800 : nn = nn - ll
3479 : END IF
3480 : END DO
3481 : END DO
3482 1659128 : ncalc = nelem - ncore
3483 : END IF
3484 :
3485 : ! readjust the occupation number of the orbitals as requested by user
3486 : ! this is done to break symmetry (bs) and bias the initial guess
3487 : ! to the pre-defined multiplicity/charge state of the atom
3488 52771 : IF (bs_occupation) THEN
3489 672 : DO is = 1, 2
3490 1216 : DO i = 1, SIZE(addel, 1)
3491 544 : ne = addel(i, is)
3492 544 : l = laddel(i, is)
3493 544 : nn = naddel(i, is) - l
3494 992 : IF (ne /= 0) THEN
3495 516 : IF (nn == 0) THEN
3496 0 : DO ii = SIZE(nelem, 2), 1, -1
3497 0 : IF (ncalc(l, ii) > 0) THEN
3498 0 : IF ((ncalc(l, ii) + ne) < 2*(2*l + 1) + 1) THEN
3499 0 : edelta(l, ii, is) = edelta(l, ii, is) + ne
3500 0 : nn = ii
3501 : ELSE
3502 0 : edelta(l, ii + 1, is) = edelta(l, ii + 1, is) + ne
3503 0 : nn = ii + 1
3504 : END IF
3505 : EXIT
3506 0 : ELSE IF (ii == 1) THEN
3507 0 : edelta(l, ii, is) = edelta(l, ii, is) + ne
3508 0 : nn = ii
3509 : END IF
3510 : END DO
3511 : ELSE
3512 516 : edelta(l, nn, is) = edelta(l, nn, is) + ne
3513 : END IF
3514 516 : IF (ncalc(l, nn) + edelta(l, nn, is) < 0) THEN
3515 0 : edelta(l, nn, is) = -ncalc(l, nn)
3516 : END IF
3517 : END IF
3518 : END DO
3519 : END DO
3520 32032 : edelta = 0.5_dp*edelta
3521 52547 : ELSE IF (magnetization /= 0.0_dp) THEN
3522 6 : dmag = 0.5_dp*ABS(magnetization)
3523 18 : DO l = 0, MIN(lmat, UBOUND(ptable(z)%e_conv, 1))
3524 18 : ll = 2*(2*l + 1)
3525 18 : ii = 0
3526 138 : DO i = 1, SIZE(ncalc, 2)
3527 126 : IF (ncalc(l, i) == 0) CYCLE
3528 24 : IF (ncalc(l, i) == ll) CYCLE
3529 24 : IF (ncalc(l, i) > dmag .AND. (ll - ncalc(l, i)) > dmag) THEN
3530 : ii = i
3531 : EXIT
3532 : END IF
3533 : END DO
3534 18 : IF (ii /= 0) THEN
3535 6 : edelta(l, ii, 1) = magnetization*0.5_dp
3536 6 : edelta(l, ii, 2) = -magnetization*0.5_dp
3537 6 : EXIT
3538 : END IF
3539 : END DO
3540 6 : IF (ii == 0) THEN
3541 : CALL cp_abort(__LOCATION__, &
3542 0 : "Magnetization value cannot be imposed for this atom type")
3543 : END IF
3544 : END IF
3545 :
3546 52771 : IF (qs_kind%ghost .OR. qs_kind%floating) THEN
3547 404 : nelem = 0
3548 404 : ncore = 0
3549 404 : ncalc = 0
3550 404 : edelta = 0.0_dp
3551 : END IF
3552 :
3553 52771 : END SUBROUTINE init_atom_electronic_state
3554 :
3555 : ! **************************************************************************************************
3556 : !> \brief ...
3557 : !> \param econf ...
3558 : !> \param z ...
3559 : !> \param ncalc ...
3560 : !> \param ncore ...
3561 : !> \param nelem ...
3562 : ! **************************************************************************************************
3563 29459 : SUBROUTINE set_pseudo_state(econf, z, ncalc, ncore, nelem)
3564 : INTEGER, DIMENSION(:), POINTER :: econf
3565 : INTEGER, INTENT(IN) :: z
3566 : INTEGER, DIMENSION(0:lmat, 10), INTENT(OUT) :: ncalc, ncore, nelem
3567 :
3568 : CHARACTER(LEN=default_string_length) :: message
3569 : INTEGER :: ii, iounit, l, ll, lmin, nc, nn
3570 : INTEGER, DIMENSION(0:lmat) :: econfx
3571 : TYPE(cp_logger_type), POINTER :: logger
3572 :
3573 29459 : NULLIFY (logger)
3574 29459 : logger => cp_get_default_logger()
3575 29459 : iounit = cp_logger_get_default_io_unit(logger)
3576 :
3577 29459 : econfx = 0
3578 83476 : econfx(0:SIZE(econf) - 1) = econf
3579 83476 : IF (SUM(econf) >= 0) THEN
3580 83408 : lmin = MIN(lmat, UBOUND(ptable(z)%e_conv, 1))
3581 : ! number of core electrons
3582 83408 : nc = z - SUM(econf)
3583 : ! setup ncore
3584 29425 : ncore = 0
3585 10840 : SELECT CASE (nc)
3586 : CASE (0)
3587 : CASE (2)
3588 10840 : ncore(0, 1) = 2
3589 : CASE (10)
3590 2448 : ncore(0, 1) = 2
3591 2448 : ncore(0, 2) = 2
3592 2448 : ncore(1, 1) = 6
3593 : CASE (18)
3594 58 : ncore(0, 1) = 2
3595 58 : ncore(0, 2) = 2
3596 58 : ncore(0, 3) = 2
3597 58 : ncore(1, 1) = 6
3598 58 : ncore(1, 2) = 6
3599 : CASE (28)
3600 22 : ncore(0, 1) = 2
3601 22 : ncore(0, 2) = 2
3602 22 : ncore(0, 3) = 2
3603 22 : ncore(1, 1) = 6
3604 22 : ncore(1, 2) = 6
3605 22 : ncore(2, 1) = 10
3606 : CASE (36)
3607 0 : ncore(0, 1) = 2
3608 0 : ncore(0, 2) = 2
3609 0 : ncore(0, 3) = 2
3610 0 : ncore(0, 4) = 2
3611 0 : ncore(1, 1) = 6
3612 0 : ncore(1, 2) = 6
3613 0 : ncore(1, 3) = 6
3614 0 : ncore(2, 1) = 10
3615 : CASE (46)
3616 72 : ncore(0, 1) = 2
3617 72 : ncore(0, 2) = 2
3618 72 : ncore(0, 3) = 2
3619 72 : ncore(0, 4) = 2
3620 72 : ncore(1, 1) = 6
3621 72 : ncore(1, 2) = 6
3622 72 : ncore(1, 3) = 6
3623 72 : ncore(2, 1) = 10
3624 72 : ncore(2, 2) = 10
3625 : CASE (54)
3626 4 : ncore(0, 1) = 2
3627 4 : ncore(0, 2) = 2
3628 4 : ncore(0, 3) = 2
3629 4 : ncore(0, 4) = 2
3630 4 : ncore(0, 5) = 2
3631 4 : ncore(1, 1) = 6
3632 4 : ncore(1, 2) = 6
3633 4 : ncore(1, 3) = 6
3634 4 : ncore(1, 4) = 6
3635 4 : ncore(2, 1) = 10
3636 4 : ncore(2, 2) = 10
3637 : CASE (60)
3638 34 : ncore(0, 1) = 2
3639 34 : ncore(0, 2) = 2
3640 34 : ncore(0, 3) = 2
3641 34 : ncore(0, 4) = 2
3642 34 : ncore(1, 1) = 6
3643 34 : ncore(1, 2) = 6
3644 34 : ncore(1, 3) = 6
3645 34 : ncore(2, 1) = 10
3646 34 : ncore(2, 2) = 10
3647 34 : ncore(3, 1) = 14
3648 : CASE (68)
3649 250 : ncore(0, 1) = 2
3650 250 : ncore(0, 2) = 2
3651 250 : ncore(0, 3) = 2
3652 250 : ncore(0, 4) = 2
3653 250 : ncore(0, 5) = 2
3654 250 : ncore(1, 1) = 6
3655 250 : ncore(1, 2) = 6
3656 250 : ncore(1, 3) = 6
3657 250 : ncore(1, 4) = 6
3658 250 : ncore(2, 1) = 10
3659 250 : ncore(2, 2) = 10
3660 250 : ncore(3, 1) = 14
3661 : CASE (78)
3662 18 : ncore(0, 1) = 2
3663 18 : ncore(0, 2) = 2
3664 18 : ncore(0, 3) = 2
3665 18 : ncore(0, 4) = 2
3666 18 : ncore(0, 5) = 2
3667 18 : ncore(1, 1) = 6
3668 18 : ncore(1, 2) = 6
3669 18 : ncore(1, 3) = 6
3670 18 : ncore(1, 4) = 6
3671 18 : ncore(2, 1) = 10
3672 18 : ncore(2, 2) = 10
3673 18 : ncore(2, 3) = 10
3674 18 : ncore(3, 1) = 14
3675 : ! 79 - 92 5f incore PP
3676 : CASE (79)
3677 0 : ncore(0, 1) = 2
3678 0 : ncore(0, 2) = 2
3679 0 : ncore(0, 3) = 2
3680 0 : ncore(0, 4) = 2
3681 0 : ncore(0, 5) = 2
3682 0 : ncore(1, 1) = 6
3683 0 : ncore(1, 2) = 6
3684 0 : ncore(1, 3) = 6
3685 0 : ncore(1, 4) = 6
3686 0 : ncore(2, 1) = 10
3687 0 : ncore(2, 2) = 10
3688 0 : ncore(2, 3) = 10
3689 0 : ncore(3, 1) = 14
3690 0 : ncore(3, 2) = 1
3691 : CASE (80)
3692 0 : ncore(0, 1) = 2
3693 0 : ncore(0, 2) = 2
3694 0 : ncore(0, 3) = 2
3695 0 : ncore(0, 4) = 2
3696 0 : ncore(0, 5) = 2
3697 0 : ncore(1, 1) = 6
3698 0 : ncore(1, 2) = 6
3699 0 : ncore(1, 3) = 6
3700 0 : ncore(1, 4) = 6
3701 0 : ncore(2, 1) = 10
3702 0 : ncore(2, 2) = 10
3703 0 : ncore(2, 3) = 10
3704 0 : ncore(3, 1) = 14
3705 0 : ncore(3, 2) = 2
3706 : CASE (81)
3707 0 : ncore(0, 1) = 2
3708 0 : ncore(0, 2) = 2
3709 0 : ncore(0, 3) = 2
3710 0 : ncore(0, 4) = 2
3711 0 : ncore(0, 5) = 2
3712 0 : ncore(1, 1) = 6
3713 0 : ncore(1, 2) = 6
3714 0 : ncore(1, 3) = 6
3715 0 : ncore(1, 4) = 6
3716 0 : ncore(2, 1) = 10
3717 0 : ncore(2, 2) = 10
3718 0 : ncore(2, 3) = 10
3719 0 : ncore(3, 1) = 14
3720 0 : ncore(3, 2) = 3
3721 : CASE (82)
3722 0 : ncore(0, 1) = 2
3723 0 : ncore(0, 2) = 2
3724 0 : ncore(0, 3) = 2
3725 0 : ncore(0, 4) = 2
3726 0 : ncore(0, 5) = 2
3727 0 : ncore(1, 1) = 6
3728 0 : ncore(1, 2) = 6
3729 0 : ncore(1, 3) = 6
3730 0 : ncore(1, 4) = 6
3731 0 : ncore(2, 1) = 10
3732 0 : ncore(2, 2) = 10
3733 0 : ncore(2, 3) = 10
3734 0 : ncore(3, 1) = 14
3735 0 : ncore(3, 2) = 4
3736 : CASE (83)
3737 0 : ncore(0, 1) = 2
3738 0 : ncore(0, 2) = 2
3739 0 : ncore(0, 3) = 2
3740 0 : ncore(0, 4) = 2
3741 0 : ncore(0, 5) = 2
3742 0 : ncore(1, 1) = 6
3743 0 : ncore(1, 2) = 6
3744 0 : ncore(1, 3) = 6
3745 0 : ncore(1, 4) = 6
3746 0 : ncore(2, 1) = 10
3747 0 : ncore(2, 2) = 10
3748 0 : ncore(2, 3) = 10
3749 0 : ncore(3, 1) = 14
3750 0 : ncore(3, 2) = 5
3751 : CASE (84)
3752 0 : ncore(0, 1) = 2
3753 0 : ncore(0, 2) = 2
3754 0 : ncore(0, 3) = 2
3755 0 : ncore(0, 4) = 2
3756 0 : ncore(0, 5) = 2
3757 0 : ncore(1, 1) = 6
3758 0 : ncore(1, 2) = 6
3759 0 : ncore(1, 3) = 6
3760 0 : ncore(1, 4) = 6
3761 0 : ncore(2, 1) = 10
3762 0 : ncore(2, 2) = 10
3763 0 : ncore(2, 3) = 10
3764 0 : ncore(3, 1) = 14
3765 0 : ncore(3, 2) = 6
3766 : CASE (85)
3767 0 : ncore(0, 1) = 2
3768 0 : ncore(0, 2) = 2
3769 0 : ncore(0, 3) = 2
3770 0 : ncore(0, 4) = 2
3771 0 : ncore(0, 5) = 2
3772 0 : ncore(1, 1) = 6
3773 0 : ncore(1, 2) = 6
3774 0 : ncore(1, 3) = 6
3775 0 : ncore(1, 4) = 6
3776 0 : ncore(2, 1) = 10
3777 0 : ncore(2, 2) = 10
3778 0 : ncore(2, 3) = 10
3779 0 : ncore(3, 1) = 14
3780 0 : ncore(3, 2) = 7
3781 : CASE (86)
3782 : ! this is not Rn core, add double assignment below
3783 0 : ncore(0, 1) = 2
3784 0 : ncore(0, 2) = 2
3785 0 : ncore(0, 3) = 2
3786 0 : ncore(0, 4) = 2
3787 0 : ncore(0, 5) = 2
3788 0 : ncore(1, 1) = 6
3789 0 : ncore(1, 2) = 6
3790 0 : ncore(1, 3) = 6
3791 0 : ncore(1, 4) = 6
3792 0 : ncore(2, 1) = 10
3793 0 : ncore(2, 2) = 10
3794 0 : ncore(2, 3) = 10
3795 0 : ncore(3, 1) = 14
3796 0 : ncore(3, 2) = 8
3797 : CASE (87)
3798 0 : ncore(0, 1) = 2
3799 0 : ncore(0, 2) = 2
3800 0 : ncore(0, 3) = 2
3801 0 : ncore(0, 4) = 2
3802 0 : ncore(0, 5) = 2
3803 0 : ncore(1, 1) = 6
3804 0 : ncore(1, 2) = 6
3805 0 : ncore(1, 3) = 6
3806 0 : ncore(1, 4) = 6
3807 0 : ncore(2, 1) = 10
3808 0 : ncore(2, 2) = 10
3809 0 : ncore(2, 3) = 10
3810 0 : ncore(3, 1) = 14
3811 0 : ncore(3, 2) = 9
3812 : CASE (88)
3813 0 : ncore(0, 1) = 2
3814 0 : ncore(0, 2) = 2
3815 0 : ncore(0, 3) = 2
3816 0 : ncore(0, 4) = 2
3817 0 : ncore(0, 5) = 2
3818 0 : ncore(1, 1) = 6
3819 0 : ncore(1, 2) = 6
3820 0 : ncore(1, 3) = 6
3821 0 : ncore(1, 4) = 6
3822 0 : ncore(2, 1) = 10
3823 0 : ncore(2, 2) = 10
3824 0 : ncore(2, 3) = 10
3825 0 : ncore(3, 1) = 14
3826 0 : ncore(3, 2) = 10
3827 : CASE (89)
3828 0 : ncore(0, 1) = 2
3829 0 : ncore(0, 2) = 2
3830 0 : ncore(0, 3) = 2
3831 0 : ncore(0, 4) = 2
3832 0 : ncore(0, 5) = 2
3833 0 : ncore(1, 1) = 6
3834 0 : ncore(1, 2) = 6
3835 0 : ncore(1, 3) = 6
3836 0 : ncore(1, 4) = 6
3837 0 : ncore(2, 1) = 10
3838 0 : ncore(2, 2) = 10
3839 0 : ncore(2, 3) = 10
3840 0 : ncore(3, 1) = 14
3841 0 : ncore(3, 2) = 11
3842 : CASE (90)
3843 0 : ncore(0, 1) = 2
3844 0 : ncore(0, 2) = 2
3845 0 : ncore(0, 3) = 2
3846 0 : ncore(0, 4) = 2
3847 0 : ncore(0, 5) = 2
3848 0 : ncore(1, 1) = 6
3849 0 : ncore(1, 2) = 6
3850 0 : ncore(1, 3) = 6
3851 0 : ncore(1, 4) = 6
3852 0 : ncore(2, 1) = 10
3853 0 : ncore(2, 2) = 10
3854 0 : ncore(2, 3) = 10
3855 0 : ncore(3, 1) = 14
3856 0 : ncore(3, 2) = 12
3857 : CASE (91)
3858 0 : ncore(0, 1) = 2
3859 0 : ncore(0, 2) = 2
3860 0 : ncore(0, 3) = 2
3861 0 : ncore(0, 4) = 2
3862 0 : ncore(0, 5) = 2
3863 0 : ncore(1, 1) = 6
3864 0 : ncore(1, 2) = 6
3865 0 : ncore(1, 3) = 6
3866 0 : ncore(1, 4) = 6
3867 0 : ncore(2, 1) = 10
3868 0 : ncore(2, 2) = 10
3869 0 : ncore(2, 3) = 10
3870 0 : ncore(3, 1) = 14
3871 0 : ncore(3, 2) = 13
3872 : CASE (92)
3873 0 : ncore(0, 1) = 2
3874 0 : ncore(0, 2) = 2
3875 0 : ncore(0, 3) = 2
3876 0 : ncore(0, 4) = 2
3877 0 : ncore(0, 5) = 2
3878 0 : ncore(1, 1) = 6
3879 0 : ncore(1, 2) = 6
3880 0 : ncore(1, 3) = 6
3881 0 : ncore(1, 4) = 6
3882 0 : ncore(2, 1) = 10
3883 0 : ncore(2, 2) = 10
3884 0 : ncore(2, 3) = 10
3885 0 : ncore(3, 1) = 14
3886 0 : ncore(3, 2) = 14
3887 : CASE DEFAULT
3888 29425 : ncore(0, 1) = -1
3889 : END SELECT
3890 : ! special cases of double assignments
3891 29425 : IF (z == 65 .AND. econfx(3) == 0) THEN
3892 : ! 4f in core for Tb
3893 4 : ncore = 0
3894 4 : ncore(0, 1) = -1
3895 : END IF
3896 : ! if there is still no core, check for special cases
3897 29425 : IF (ncore(0, 1) <= 0) THEN
3898 15683 : IF (z >= 58 .AND. z <= 71) THEN
3899 : ! 4f-in-core PPs for lanthanides
3900 280 : nc = z - SUM(econf)
3901 : ! setup ncore
3902 56 : ncore = 0
3903 0 : SELECT CASE (nc)
3904 : CASE (29:42)
3905 0 : ncore(0, 1) = 2
3906 0 : ncore(0, 2) = 2
3907 0 : ncore(0, 3) = 2
3908 0 : ncore(1, 1) = 6
3909 0 : ncore(1, 2) = 6
3910 0 : ncore(2, 1) = 10
3911 0 : ncore(3, 1) = nc - 28
3912 : message = "A small-core pseudopotential with 4f-in-core is used for the lanthanide "// &
3913 0 : TRIM(ptable(z)%symbol)
3914 0 : CPHINT(TRIM(message))
3915 : CASE (47:60)
3916 56 : ncore(0, 1) = 2
3917 56 : ncore(0, 2) = 2
3918 56 : ncore(0, 3) = 2
3919 56 : ncore(0, 4) = 2
3920 56 : ncore(1, 1) = 6
3921 56 : ncore(1, 2) = 6
3922 56 : ncore(1, 3) = 6
3923 56 : ncore(2, 1) = 10
3924 56 : ncore(2, 2) = 10
3925 56 : ncore(3, 1) = nc - 46
3926 : message = "A medium-core pseudopotential with 4f-in-core is used for the lanthanide "// &
3927 56 : TRIM(ptable(z)%symbol)
3928 56 : CPHINT(TRIM(message))
3929 : CASE DEFAULT
3930 56 : ncore(0, 1) = -1
3931 : END SELECT
3932 : END IF
3933 : END IF
3934 : ! if the core is established, finish the setup
3935 29425 : IF (ncore(0, 1) >= 0) THEN
3936 147125 : DO l = 0, lmin
3937 117700 : ll = 2*(2*l + 1)
3938 1294700 : nn = SUM(ncore(l, :)) + econfx(l)
3939 117700 : ii = 0
3940 29425 : DO
3941 140594 : ii = ii + 1
3942 140594 : IF (nn <= ll) THEN
3943 117700 : nelem(l, ii) = nn
3944 : EXIT
3945 : ELSE
3946 22894 : nelem(l, ii) = ll
3947 22894 : nn = nn - ll
3948 : END IF
3949 : END DO
3950 : END DO
3951 2089175 : ncalc = nelem - ncore
3952 : ELSE
3953 : ! test for compatibility of valence occupation and full atomic occupation
3954 0 : IF (iounit > 0) THEN
3955 0 : WRITE (iounit, "(/,A,A2)") "WARNING: Core states irregular for atom type ", ptable(z)%symbol
3956 0 : WRITE (iounit, "(A,10I3)") "WARNING: Redefine ELEC_CONF in the KIND section"
3957 0 : CPABORT("Incompatible Atomic Occupations Detected")
3958 : END IF
3959 : END IF
3960 : ELSE
3961 34 : lmin = MIN(lmat, UBOUND(ptable(z)%e_conv, 1))
3962 34 : ncore = 0
3963 34 : ncalc = 0
3964 170 : DO l = 0, lmin
3965 136 : ll = 2*(2*l + 1)
3966 136 : nn = ABS(econfx(l))
3967 136 : ii = 0
3968 34 : DO
3969 136 : ii = ii + 1
3970 136 : IF (nn <= ll) THEN
3971 136 : ncalc(l, ii) = -nn
3972 : EXIT
3973 : ELSE
3974 0 : ncalc(l, ii) = -ll
3975 0 : nn = nn - ll
3976 : END IF
3977 : END DO
3978 : END DO
3979 34 : nelem = ncalc
3980 : END IF
3981 :
3982 29459 : END SUBROUTINE set_pseudo_state
3983 :
3984 : ! **************************************************************************************************
3985 : !> \brief finds if a given qs run needs to use nlcc
3986 : !> \param qs_kind_set ...
3987 : !> \return ...
3988 : ! **************************************************************************************************
3989 44380 : FUNCTION has_nlcc(qs_kind_set) RESULT(nlcc)
3990 :
3991 : TYPE(qs_kind_type), DIMENSION(:) :: qs_kind_set
3992 : LOGICAL :: nlcc
3993 :
3994 : INTEGER :: ikind
3995 : LOGICAL :: nlcc_present
3996 : TYPE(gth_potential_type), POINTER :: gth_potential
3997 : TYPE(sgp_potential_type), POINTER :: sgp_potential
3998 :
3999 44380 : nlcc = .FALSE.
4000 :
4001 129911 : DO ikind = 1, SIZE(qs_kind_set)
4002 85531 : CALL get_qs_kind(qs_kind_set(ikind), gth_potential=gth_potential, sgp_potential=sgp_potential)
4003 129911 : IF (ASSOCIATED(gth_potential)) THEN
4004 53927 : CALL get_potential(potential=gth_potential, nlcc_present=nlcc_present)
4005 53927 : nlcc = nlcc .OR. nlcc_present
4006 31604 : ELSE IF (ASSOCIATED(sgp_potential)) THEN
4007 836 : CALL get_potential(potential=sgp_potential, has_nlcc=nlcc_present)
4008 836 : nlcc = nlcc .OR. nlcc_present
4009 : END IF
4010 : END DO
4011 :
4012 44380 : END FUNCTION has_nlcc
4013 :
4014 : ! **************************************************************************************************
4015 :
4016 0 : END MODULE qs_kind_types
|
