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