Line data Source code
1 : !--------------------------------------------------------------------------------------------------!
2 : ! CP2K: A general program to perform molecular dynamics simulations !
3 : ! Copyright 2000-2026 CP2K developers group <https://cp2k.org> !
4 : ! !
5 : ! SPDX-License-Identifier: GPL-2.0-or-later !
6 : !--------------------------------------------------------------------------------------------------!
7 :
8 : ! **************************************************************************************************
9 : !> \brief Utilities to set up the control types
10 : ! **************************************************************************************************
11 : MODULE cp_control_utils
12 : USE bibliography, ONLY: &
13 : Andreussi2012, Andreussi2019, Chai2025a, Dewar1977, Dewar1985, Elstner1998, Fattebert2002, &
14 : Grimme2017, Hu2007, Katbashev2025, Krack2000, Lippert1997, Lippert1999, Porezag1995, &
15 : Pracht2019, Repasky2002, Rocha2006, Schenter2008, Seifert1996, Souza2002, Stengel2009, &
16 : Stewart1989, Stewart2007, Thiel1992, Umari2002, VanVoorhis2015, VandeVondele2005a, &
17 : VandeVondele2005b, Yin2017, Zhechkov2005, cite_reference
18 : USE cell_types, ONLY: cell_transform_input_cartesian,&
19 : cell_type
20 : USE cp_control_types, ONLY: &
21 : admm_control_create, admm_control_type, ddapc_control_create, ddapc_restraint_type, &
22 : dft_control_create, dft_control_type, efield_type, expot_control_create, &
23 : maxwell_control_create, qs_control_type, rixs_control_type, tddfpt2_control_type, &
24 : xtb_control_type, xtb_reference_cli_type
25 : USE cp_files, ONLY: close_file,&
26 : open_file
27 : USE cp_log_handling, ONLY: cp_get_default_logger,&
28 : cp_logger_type
29 : USE cp_output_handling, ONLY: cp_print_key_finished_output,&
30 : cp_print_key_unit_nr
31 : USE cp_parser_methods, ONLY: parser_read_line
32 : USE cp_parser_types, ONLY: cp_parser_type,&
33 : parser_create,&
34 : parser_release,&
35 : parser_reset
36 : USE cp_spline_utils, ONLY: pw_interp
37 : USE cp_units, ONLY: cp_unit_from_cp2k,&
38 : cp_unit_to_cp2k
39 : USE eeq_input, ONLY: read_eeq_param
40 : USE force_fields_input, ONLY: read_gp_section
41 : USE input_constants, ONLY: &
42 : admm1_type, admm2_type, admmp_type, admmq_type, admms_type, constant_env, custom_env, &
43 : do_admm_aux_exch_func_bee, do_admm_aux_exch_func_bee_libxc, do_admm_aux_exch_func_default, &
44 : do_admm_aux_exch_func_default_libxc, do_admm_aux_exch_func_none, &
45 : do_admm_aux_exch_func_opt, do_admm_aux_exch_func_opt_libxc, do_admm_aux_exch_func_pbex, &
46 : do_admm_aux_exch_func_pbex_libxc, do_admm_aux_exch_func_sx_libxc, &
47 : do_admm_basis_projection, do_admm_blocked_projection, do_admm_blocking_purify_full, &
48 : do_admm_charge_constrained_projection, do_admm_exch_scaling_merlot, &
49 : do_admm_exch_scaling_none, do_admm_purify_cauchy, do_admm_purify_cauchy_subspace, &
50 : do_admm_purify_mcweeny, do_admm_purify_mo_diag, do_admm_purify_mo_no_diag, &
51 : do_admm_purify_none, do_admm_purify_none_dm, do_ddapc_constraint, do_ddapc_restraint, &
52 : do_method_am1, do_method_dftb, do_method_gapw, do_method_gapw_xc, do_method_gpw, &
53 : do_method_lrigpw, do_method_mndo, do_method_mndod, do_method_ofgpw, do_method_pdg, &
54 : do_method_pm3, do_method_pm6, do_method_pm6fm, do_method_pnnl, do_method_rigpw, &
55 : do_method_rm1, do_method_xtb, do_pwgrid_ns_fullspace, do_pwgrid_ns_halfspace, &
56 : do_pwgrid_spherical, do_s2_constraint, do_s2_restraint, do_se_is_kdso, do_se_is_kdso_d, &
57 : do_se_is_slater, do_se_lr_ewald, do_se_lr_ewald_gks, do_se_lr_ewald_r3, do_se_lr_none, &
58 : gapw_1c_large, gapw_1c_medium, gapw_1c_orb, gapw_1c_small, gapw_1c_very_large, &
59 : gaussian_env, gfn1xtb, gfn_tblite, kg_tnadd_embed, kg_tnadd_embed_ri, no_admm_type, &
60 : numerical, ramp_env, real_time_propagation, rtp_method_bse, rtp_method_bse_linearized, &
61 : sccs_andreussi, sccs_derivative_cd3, sccs_derivative_cd5, sccs_derivative_cd7, &
62 : sccs_derivative_fft, sccs_fattebert_gygi, sccs_saa_andreussi, sic_ad, sic_eo, &
63 : sic_list_all, sic_list_unpaired, sic_mauri_spz, sic_mauri_us, sic_none, slater, &
64 : tblite_cli_born_kernel_auto, tblite_cli_solution_state_gsolv, tblite_cli_solvation_alpb, &
65 : tblite_cli_solvation_cpcm, tblite_cli_solvation_gb, tblite_cli_solvation_gbe, &
66 : tblite_cli_solvation_gbsa, tblite_guess_ceh, tblite_mixer_memory_inherit, &
67 : tblite_scc_mixer_auto, tblite_scc_mixer_cp2k, tblite_scc_mixer_none, &
68 : tblite_scc_mixer_tblite, tblite_solver_gvd, tblite_solver_gvr, tddfpt_dipole_length, &
69 : tddfpt_kernel_stda, use_mom_ref_user, xtb_vdw_type_d3, xtb_vdw_type_d4, xtb_vdw_type_none
70 : USE input_cp2k_check, ONLY: xc_functionals_expand
71 : USE input_cp2k_dft, ONLY: create_dft_section
72 : USE input_enumeration_types, ONLY: enum_i2c,&
73 : enumeration_type
74 : USE input_keyword_types, ONLY: keyword_get,&
75 : keyword_type
76 : USE input_section_types, ONLY: &
77 : section_get_ival, section_get_keyword, section_release, section_type, section_vals_get, &
78 : section_vals_get_subs_vals, section_vals_type, section_vals_val_get, section_vals_val_set
79 : USE kinds, ONLY: default_path_length,&
80 : default_string_length,&
81 : dp
82 : USE mathconstants, ONLY: fourpi
83 : USE pair_potential_types, ONLY: pair_potential_reallocate
84 : USE periodic_table, ONLY: get_ptable_info
85 : USE qs_cdft_utils, ONLY: read_cdft_control_section
86 : USE smeagol_control_types, ONLY: read_smeagol_control
87 : USE string_utilities, ONLY: uppercase
88 : USE util, ONLY: sort
89 : USE xas_tdp_types, ONLY: read_xas_tdp_control
90 : USE xc, ONLY: xc_uses_kinetic_energy_density,&
91 : xc_uses_norm_drho
92 : USE xc_input_constants, ONLY: xc_deriv_collocate
93 : USE xc_write_output, ONLY: xc_write
94 : #include "./base/base_uses.f90"
95 :
96 : IMPLICIT NONE
97 :
98 : PRIVATE
99 :
100 : CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'cp_control_utils'
101 :
102 : PUBLIC :: read_dft_control, &
103 : read_rixs_control, &
104 : read_mgrid_section, &
105 : read_qs_section, &
106 : read_tddfpt2_control, &
107 : write_dft_control, &
108 : write_qs_control, &
109 : write_admm_control, &
110 : read_ddapc_section
111 : CONTAINS
112 :
113 : ! **************************************************************************************************
114 : !> \brief ...
115 : !> \param dft_control ...
116 : !> \param dft_section ...
117 : !> \param cell ...
118 : ! **************************************************************************************************
119 164556 : SUBROUTINE read_dft_control(dft_control, dft_section, cell)
120 : TYPE(dft_control_type), POINTER :: dft_control
121 : TYPE(section_vals_type), POINTER :: dft_section
122 : TYPE(cell_type), OPTIONAL, POINTER :: cell
123 :
124 : CHARACTER(len=default_path_length) :: basis_set_file_name, gauxc_model_name, &
125 : intensities_file_name, &
126 : potential_file_name
127 : CHARACTER(LEN=default_string_length), &
128 9142 : DIMENSION(:), POINTER :: tmpstringlist
129 : INTEGER :: admmtype, irep, isize, kg_tnadd_method, &
130 : method_id, nrep, xc_deriv_method_id
131 : LOGICAL :: at_end, do_hfx, do_ot, do_rpa_admm, do_rtp, exopt1, exopt2, exopt3, explicit, &
132 : is_present, l_param, local_moment_possible, native_skala_grid, not_SE, was_present
133 : REAL(KIND=dp) :: density_cut, gradient_cut, tau_cut
134 9142 : REAL(KIND=dp), DIMENSION(:), POINTER :: pol
135 : TYPE(cp_logger_type), POINTER :: logger
136 : TYPE(cp_parser_type) :: parser
137 : TYPE(section_vals_type), POINTER :: hairy_probes_section, hfx_section, kg_xc_fun_section, &
138 : kg_xc_section, maxwell_section, sccs_section, scf_section, tmp_section, xc_fun_section, &
139 : xc_gauxc_subsection, xc_section
140 :
141 9142 : was_present = .FALSE.
142 :
143 9142 : logger => cp_get_default_logger()
144 :
145 9142 : NULLIFY (kg_xc_fun_section, kg_xc_section, tmp_section, xc_fun_section, xc_section, xc_gauxc_subsection)
146 45710 : ALLOCATE (dft_control)
147 9142 : CALL dft_control_create(dft_control)
148 : ! determine wheather this is a semiempirical or DFTB run
149 : ! --> (no XC section needs to be provided)
150 9142 : not_SE = .TRUE.
151 9142 : CALL section_vals_val_get(dft_section, "QS%METHOD", i_val=method_id)
152 2534 : SELECT CASE (method_id)
153 : CASE (do_method_dftb, do_method_xtb, do_method_mndo, do_method_am1, do_method_pm3, do_method_pnnl, &
154 : do_method_pm6, do_method_pm6fm, do_method_pdg, do_method_rm1, do_method_mndod)
155 9142 : not_SE = .FALSE.
156 : END SELECT
157 : ! Check for XC section and XC_FUNCTIONAL section
158 9142 : xc_section => section_vals_get_subs_vals(dft_section, "XC")
159 9142 : CALL section_vals_get(xc_section, explicit=is_present)
160 9142 : IF (.NOT. is_present .AND. not_SE) THEN
161 0 : CPABORT("XC section missing.")
162 : END IF
163 9142 : IF (is_present) THEN
164 6624 : CALL section_vals_val_get(xc_section, "density_cutoff", r_val=density_cut)
165 6624 : CALL section_vals_val_get(xc_section, "gradient_cutoff", r_val=gradient_cut)
166 6624 : CALL section_vals_val_get(xc_section, "tau_cutoff", r_val=tau_cut)
167 : ! Perform numerical stability checks and possibly correct the issues
168 6624 : IF (density_cut <= EPSILON(0.0_dp)*100.0_dp) THEN
169 : CALL cp_warn(__LOCATION__, &
170 : "DENSITY_CUTOFF lower than 100*EPSILON, where EPSILON is the machine precision. "// &
171 0 : "This may lead to numerical problems. Setting up shake_tol to 100*EPSILON! ")
172 : END IF
173 6624 : density_cut = MAX(EPSILON(0.0_dp)*100.0_dp, density_cut)
174 6624 : IF (gradient_cut <= EPSILON(0.0_dp)*100.0_dp) THEN
175 : CALL cp_warn(__LOCATION__, &
176 : "GRADIENT_CUTOFF lower than 100*EPSILON, where EPSILON is the machine precision. "// &
177 0 : "This may lead to numerical problems. Setting up shake_tol to 100*EPSILON! ")
178 : END IF
179 6624 : gradient_cut = MAX(EPSILON(0.0_dp)*100.0_dp, gradient_cut)
180 6624 : IF (tau_cut <= EPSILON(0.0_dp)*100.0_dp) THEN
181 : CALL cp_warn(__LOCATION__, &
182 : "TAU_CUTOFF lower than 100*EPSILON, where EPSILON is the machine precision. "// &
183 0 : "This may lead to numerical problems. Setting up shake_tol to 100*EPSILON! ")
184 : END IF
185 6624 : tau_cut = MAX(EPSILON(0.0_dp)*100.0_dp, tau_cut)
186 6624 : CALL section_vals_val_set(xc_section, "density_cutoff", r_val=density_cut)
187 6624 : CALL section_vals_val_set(xc_section, "gradient_cutoff", r_val=gradient_cut)
188 6624 : CALL section_vals_val_set(xc_section, "tau_cutoff", r_val=tau_cut)
189 : END IF
190 9142 : xc_fun_section => section_vals_get_subs_vals(xc_section, "XC_FUNCTIONAL")
191 9142 : CALL section_vals_get(xc_fun_section, explicit=is_present)
192 9142 : IF (.NOT. is_present .AND. not_SE) THEN
193 0 : CPABORT("XC_FUNCTIONAL section missing.")
194 : END IF
195 :
196 9142 : dft_control%use_gauxc = .FALSE.
197 9142 : IF (is_present) THEN
198 6624 : xc_gauxc_subsection => section_vals_get_subs_vals(xc_fun_section, "GAUXC")
199 6624 : CALL section_vals_get(xc_gauxc_subsection, explicit=dft_control%use_gauxc)
200 : END IF
201 :
202 9142 : scf_section => section_vals_get_subs_vals(dft_section, "SCF")
203 9142 : CALL section_vals_val_get(dft_section, "UKS", l_val=dft_control%uks)
204 9142 : CALL section_vals_val_get(dft_section, "ROKS", l_val=dft_control%roks)
205 9142 : IF (dft_control%uks .OR. dft_control%roks) THEN
206 1937 : dft_control%nspins = 2
207 : ELSE
208 7205 : dft_control%nspins = 1
209 : END IF
210 :
211 9142 : dft_control%lsd = (dft_control%nspins > 1)
212 9142 : dft_control%use_kinetic_energy_density = xc_uses_kinetic_energy_density(xc_fun_section, dft_control%lsd)
213 9142 : IF (dft_control%use_gauxc) THEN
214 : native_skala_grid = .FALSE.
215 142 : CALL section_vals_val_get(xc_gauxc_subsection, "NATIVE_GRID", l_val=native_skala_grid)
216 142 : CALL section_vals_val_get(xc_gauxc_subsection, "MODEL", c_val=gauxc_model_name)
217 142 : gauxc_model_name = ADJUSTL(gauxc_model_name)
218 142 : CALL uppercase(gauxc_model_name)
219 142 : IF (native_skala_grid .OR. &
220 : (TRIM(gauxc_model_name) /= "" .AND. TRIM(gauxc_model_name) /= "NONE")) THEN
221 134 : dft_control%use_kinetic_energy_density = .TRUE.
222 : END IF
223 : END IF
224 9142 : tmp_section => section_vals_get_subs_vals(dft_section, "KG_METHOD")
225 9142 : CALL section_vals_get(tmp_section, explicit=explicit)
226 9142 : IF (explicit) THEN
227 82 : CALL section_vals_val_get(tmp_section, "TNADD_METHOD", i_val=kg_tnadd_method)
228 82 : IF (kg_tnadd_method == kg_tnadd_embed .OR. kg_tnadd_method == kg_tnadd_embed_ri) THEN
229 64 : kg_xc_section => section_vals_get_subs_vals(tmp_section, "XC")
230 64 : kg_xc_fun_section => section_vals_get_subs_vals(kg_xc_section, "XC_FUNCTIONAL")
231 64 : CALL section_vals_get(kg_xc_fun_section, explicit=is_present)
232 64 : IF (is_present) THEN
233 : dft_control%use_kinetic_energy_density = dft_control%use_kinetic_energy_density .OR. &
234 : xc_uses_kinetic_energy_density(kg_xc_fun_section, &
235 68 : dft_control%lsd)
236 : END IF
237 : END IF
238 : END IF
239 :
240 9142 : xc_deriv_method_id = section_get_ival(xc_section, "XC_GRID%XC_DERIV")
241 : dft_control%drho_by_collocation = (xc_uses_norm_drho(xc_fun_section, dft_control%lsd) &
242 9142 : .AND. (xc_deriv_method_id == xc_deriv_collocate))
243 9142 : IF (dft_control%drho_by_collocation) THEN
244 0 : CPABORT("derivatives by collocation not implemented")
245 : END IF
246 :
247 : ! Automatic auxiliary basis set generation
248 9142 : CALL section_vals_val_get(dft_section, "AUTO_BASIS", n_rep_val=nrep)
249 18284 : DO irep = 1, nrep
250 9142 : CALL section_vals_val_get(dft_section, "AUTO_BASIS", i_rep_val=irep, c_vals=tmpstringlist)
251 18284 : IF (SIZE(tmpstringlist) == 2) THEN
252 9142 : CALL uppercase(tmpstringlist(2))
253 18128 : SELECT CASE (tmpstringlist(2))
254 : CASE ("X")
255 8986 : SELECT CASE (tmpstringlist(1))
256 : CASE ("X")
257 : ! Do nothing
258 : CASE DEFAULT
259 : CALL cp_abort(__LOCATION__, &
260 : "AUTO_BASIS: the size <X> is invalid for the "// &
261 : "type <"//TRIM(ADJUSTL(tmpstringlist(1)))//">; "// &
262 : "use one of SMALL, MEDIUM, LARGE, HUGE for "// &
263 : "the size. The syntax AUTO_BASIS X X is a "// &
264 : "reserved case for using NO automatically "// &
265 8986 : "generated basis sets.")
266 : END SELECT
267 : CASE ("SMALL")
268 54 : isize = 0
269 : CASE ("MEDIUM")
270 54 : isize = 1
271 : CASE ("LARGE")
272 0 : isize = 2
273 : CASE ("HUGE")
274 8 : isize = 3
275 : CASE DEFAULT
276 9142 : CPWARN("Unknown basis size in AUTO_BASIS keyword:"//TRIM(tmpstringlist(1)))
277 : END SELECT
278 : !
279 9144 : SELECT CASE (tmpstringlist(1))
280 : CASE ("X")
281 : CASE ("RI_AUX")
282 2 : dft_control%auto_basis_ri_aux = isize
283 : CASE ("AUX_FIT")
284 0 : dft_control%auto_basis_aux_fit = isize
285 : CASE ("LRI_AUX")
286 0 : dft_control%auto_basis_lri_aux = isize
287 : CASE ("P_LRI_AUX")
288 0 : dft_control%auto_basis_p_lri_aux = isize
289 : CASE ("RI_HXC")
290 0 : dft_control%auto_basis_ri_hxc = isize
291 : CASE ("RI_XAS")
292 64 : dft_control%auto_basis_ri_xas = isize
293 : CASE ("RI_HFX")
294 90 : dft_control%auto_basis_ri_hfx = isize
295 : CASE DEFAULT
296 9142 : CPWARN("Unknown basis type in AUTO_BASIS keyword:"//TRIM(tmpstringlist(1)))
297 : END SELECT
298 : ELSE
299 : CALL cp_abort(__LOCATION__, &
300 0 : "AUTO_BASIS keyword in &DFT section has a wrong number of arguments.")
301 : END IF
302 : END DO
303 :
304 : !! check if we do wavefunction fitting
305 9142 : tmp_section => section_vals_get_subs_vals(dft_section, "AUXILIARY_DENSITY_MATRIX_METHOD")
306 9142 : CALL section_vals_get(tmp_section, explicit=is_present)
307 : !
308 9142 : hfx_section => section_vals_get_subs_vals(xc_section, "HF")
309 9142 : CALL section_vals_get(hfx_section, explicit=do_hfx)
310 9142 : CALL section_vals_val_get(xc_section, "WF_CORRELATION%RI_RPA%ADMM", l_val=do_rpa_admm)
311 9142 : is_present = is_present .AND. (do_hfx .OR. do_rpa_admm)
312 : !
313 9142 : dft_control%do_admm = is_present
314 9142 : dft_control%do_admm_mo = .FALSE.
315 9142 : dft_control%do_admm_dm = .FALSE.
316 9142 : IF (is_present) THEN
317 : do_ot = .FALSE.
318 524 : CALL section_vals_val_get(scf_section, "OT%_SECTION_PARAMETERS_", l_val=do_ot)
319 524 : CALL admm_control_create(dft_control%admm_control)
320 :
321 524 : CALL section_vals_val_get(dft_section, "AUXILIARY_DENSITY_MATRIX_METHOD%ADMM_TYPE", i_val=admmtype)
322 524 : CALL section_vals_val_get(dft_section, "AUXILIARY_DENSITY_MATRIX_METHOD%ADMM_PURIFICATION_METHOD", explicit=exopt1)
323 524 : CALL section_vals_val_get(dft_section, "AUXILIARY_DENSITY_MATRIX_METHOD%METHOD", explicit=exopt2)
324 524 : CALL section_vals_val_get(dft_section, "AUXILIARY_DENSITY_MATRIX_METHOD%EXCH_SCALING_MODEL", explicit=exopt3)
325 524 : dft_control%admm_control%admm_type = admmtype
326 506 : SELECT CASE (admmtype)
327 : CASE (no_admm_type)
328 506 : CALL section_vals_val_get(dft_section, "AUXILIARY_DENSITY_MATRIX_METHOD%ADMM_PURIFICATION_METHOD", i_val=method_id)
329 506 : dft_control%admm_control%purification_method = method_id
330 506 : CALL section_vals_val_get(dft_section, "AUXILIARY_DENSITY_MATRIX_METHOD%METHOD", i_val=method_id)
331 506 : dft_control%admm_control%method = method_id
332 506 : CALL section_vals_val_get(dft_section, "AUXILIARY_DENSITY_MATRIX_METHOD%EXCH_SCALING_MODEL", i_val=method_id)
333 506 : dft_control%admm_control%scaling_model = method_id
334 : CASE (admm1_type)
335 : ! METHOD BASIS_PROJECTION
336 : ! ADMM_PURIFICATION_METHOD choose
337 : ! EXCH_SCALING_MODEL NONE
338 4 : CALL section_vals_val_get(dft_section, "AUXILIARY_DENSITY_MATRIX_METHOD%ADMM_PURIFICATION_METHOD", i_val=method_id)
339 4 : dft_control%admm_control%purification_method = method_id
340 4 : dft_control%admm_control%method = do_admm_basis_projection
341 4 : dft_control%admm_control%scaling_model = do_admm_exch_scaling_none
342 : CASE (admm2_type)
343 : ! METHOD BASIS_PROJECTION
344 : ! ADMM_PURIFICATION_METHOD NONE
345 : ! EXCH_SCALING_MODEL NONE
346 2 : dft_control%admm_control%purification_method = do_admm_purify_none
347 2 : dft_control%admm_control%method = do_admm_basis_projection
348 2 : dft_control%admm_control%scaling_model = do_admm_exch_scaling_none
349 : CASE (admms_type)
350 : ! ADMM_PURIFICATION_METHOD NONE
351 : ! METHOD CHARGE_CONSTRAINED_PROJECTION
352 : ! EXCH_SCALING_MODEL MERLOT
353 8 : dft_control%admm_control%purification_method = do_admm_purify_none
354 8 : dft_control%admm_control%method = do_admm_charge_constrained_projection
355 8 : dft_control%admm_control%scaling_model = do_admm_exch_scaling_merlot
356 : CASE (admmp_type)
357 : ! ADMM_PURIFICATION_METHOD NONE
358 : ! METHOD BASIS_PROJECTION
359 : ! EXCH_SCALING_MODEL MERLOT
360 2 : dft_control%admm_control%purification_method = do_admm_purify_none
361 2 : dft_control%admm_control%method = do_admm_basis_projection
362 2 : dft_control%admm_control%scaling_model = do_admm_exch_scaling_merlot
363 : CASE (admmq_type)
364 : ! ADMM_PURIFICATION_METHOD NONE
365 : ! METHOD CHARGE_CONSTRAINED_PROJECTION
366 : ! EXCH_SCALING_MODEL NONE
367 2 : dft_control%admm_control%purification_method = do_admm_purify_none
368 2 : dft_control%admm_control%method = do_admm_charge_constrained_projection
369 2 : dft_control%admm_control%scaling_model = do_admm_exch_scaling_none
370 : CASE DEFAULT
371 : CALL cp_abort(__LOCATION__, &
372 524 : "ADMM_TYPE keyword in &AUXILIARY_DENSITY_MATRIX_METHOD section has a wrong value.")
373 : END SELECT
374 :
375 : CALL section_vals_val_get(dft_section, "AUXILIARY_DENSITY_MATRIX_METHOD%EPS_FILTER", &
376 524 : r_val=dft_control%admm_control%eps_filter)
377 :
378 524 : CALL section_vals_val_get(dft_section, "AUXILIARY_DENSITY_MATRIX_METHOD%EXCH_CORRECTION_FUNC", i_val=method_id)
379 524 : dft_control%admm_control%aux_exch_func = method_id
380 :
381 : ! parameters for X functional
382 524 : dft_control%admm_control%aux_exch_func_param = .FALSE.
383 : CALL section_vals_val_get(dft_section, "AUXILIARY_DENSITY_MATRIX_METHOD%OPTX_A1", explicit=explicit, &
384 524 : r_val=dft_control%admm_control%aux_x_param(1))
385 524 : IF (explicit) dft_control%admm_control%aux_exch_func_param = .TRUE.
386 : CALL section_vals_val_get(dft_section, "AUXILIARY_DENSITY_MATRIX_METHOD%OPTX_A2", explicit=explicit, &
387 524 : r_val=dft_control%admm_control%aux_x_param(2))
388 524 : IF (explicit) dft_control%admm_control%aux_exch_func_param = .TRUE.
389 : CALL section_vals_val_get(dft_section, "AUXILIARY_DENSITY_MATRIX_METHOD%OPTX_GAMMA", explicit=explicit, &
390 524 : r_val=dft_control%admm_control%aux_x_param(3))
391 524 : IF (explicit) dft_control%admm_control%aux_exch_func_param = .TRUE.
392 :
393 524 : CALL read_admm_block_list(dft_control%admm_control, dft_section)
394 :
395 : ! check for double assignments
396 2 : SELECT CASE (admmtype)
397 : CASE (admm2_type)
398 2 : IF (exopt2) CALL cp_warn(__LOCATION__, &
399 0 : "Value of ADMM_PURIFICATION_METHOD keyword will be overwritten with ADMM_TYPE selections.")
400 2 : IF (exopt3) CALL cp_warn(__LOCATION__, &
401 0 : "Value of EXCH_SCALING_MODEL keyword will be overwritten with ADMM_TYPE selections.")
402 : CASE (admm1_type, admms_type, admmp_type, admmq_type)
403 16 : IF (exopt1) CALL cp_warn(__LOCATION__, &
404 2 : "Value of METHOD keyword will be overwritten with ADMM_TYPE selections.")
405 16 : IF (exopt2) CALL cp_warn(__LOCATION__, &
406 2 : "Value of METHOD keyword will be overwritten with ADMM_TYPE selections.")
407 16 : IF (exopt3) CALL cp_warn(__LOCATION__, &
408 524 : "Value of EXCH_SCALING_MODEL keyword will be overwritten with ADMM_TYPE selections.")
409 : END SELECT
410 :
411 : ! In the case of charge-constrained projection (e.g. according to Merlot),
412 : ! there is no purification needed and hence, do_admm_purify_none has to be set.
413 :
414 : IF ((dft_control%admm_control%method == do_admm_blocking_purify_full .OR. &
415 : dft_control%admm_control%method == do_admm_blocked_projection) &
416 524 : .AND. dft_control%admm_control%scaling_model == do_admm_exch_scaling_merlot) THEN
417 0 : CPABORT("ADMM: Blocking and Merlot scaling are mutually exclusive.")
418 : END IF
419 :
420 524 : IF (dft_control%admm_control%method == do_admm_charge_constrained_projection .AND. &
421 : dft_control%admm_control%purification_method /= do_admm_purify_none) THEN
422 : CALL cp_abort(__LOCATION__, &
423 : "ADMM: In the case of METHOD=CHARGE_CONSTRAINED_PROJECTION, "// &
424 0 : "ADMM_PURIFICATION_METHOD=NONE has to be set.")
425 : END IF
426 :
427 524 : IF (dft_control%admm_control%purification_method == do_admm_purify_mo_diag .OR. &
428 : dft_control%admm_control%purification_method == do_admm_purify_mo_no_diag) THEN
429 62 : IF (dft_control%admm_control%method /= do_admm_basis_projection) THEN
430 0 : CPABORT("ADMM: Chosen purification requires BASIS_PROJECTION")
431 : END IF
432 :
433 62 : IF (.NOT. do_ot) CPABORT("ADMM: MO-based purification requires OT.")
434 : END IF
435 :
436 524 : IF (dft_control%admm_control%purification_method == do_admm_purify_none_dm .OR. &
437 : dft_control%admm_control%purification_method == do_admm_purify_mcweeny) THEN
438 14 : dft_control%do_admm_dm = .TRUE.
439 : ELSE
440 510 : dft_control%do_admm_mo = .TRUE.
441 : END IF
442 : END IF
443 :
444 : ! Set restricted to true, if both OT and ROKS are requested
445 : !MK in principle dft_control%restricted could be dropped completely like the
446 : !MK input key by using only dft_control%roks now
447 9142 : CALL section_vals_val_get(scf_section, "OT%_SECTION_PARAMETERS_", l_val=l_param)
448 9142 : dft_control%restricted = (dft_control%roks .AND. l_param)
449 :
450 9142 : CALL section_vals_val_get(dft_section, "CHARGE", i_val=dft_control%charge)
451 9142 : CALL section_vals_val_get(dft_section, "MULTIPLICITY", i_val=dft_control%multiplicity)
452 9142 : CALL section_vals_val_get(dft_section, "RELAX_MULTIPLICITY", r_val=dft_control%relax_multiplicity)
453 9142 : IF (dft_control%relax_multiplicity > 0.0_dp) THEN
454 10 : IF (.NOT. dft_control%uks) THEN
455 : CALL cp_abort(__LOCATION__, "The option RELAX_MULTIPLICITY is only valid for "// &
456 0 : "unrestricted Kohn-Sham (UKS) calculations")
457 : END IF
458 : END IF
459 :
460 : !Read the HAIR PROBES input section if present
461 9142 : hairy_probes_section => section_vals_get_subs_vals(dft_section, "HAIRY_PROBES")
462 9142 : CALL section_vals_get(hairy_probes_section, n_repetition=nrep, explicit=is_present)
463 :
464 9142 : IF (is_present) THEN
465 4 : dft_control%hairy_probes = .TRUE.
466 20 : ALLOCATE (dft_control%probe(nrep))
467 4 : CALL read_hairy_probes_sections(dft_control, hairy_probes_section)
468 : END IF
469 :
470 : ! check for the presence of the low spin roks section
471 9142 : tmp_section => section_vals_get_subs_vals(dft_section, "LOW_SPIN_ROKS")
472 9142 : CALL section_vals_get(tmp_section, explicit=dft_control%low_spin_roks)
473 :
474 9142 : dft_control%sic_method_id = sic_none
475 9142 : dft_control%sic_scaling_a = 1.0_dp
476 9142 : dft_control%sic_scaling_b = 1.0_dp
477 :
478 : ! DFT+U
479 9142 : dft_control%dft_plus_u = .FALSE.
480 9142 : CALL section_vals_val_get(dft_section, "PLUS_U_METHOD", i_val=method_id)
481 9142 : dft_control%plus_u_method_id = method_id
482 :
483 : ! Minimum tracking linear response U and J
484 9142 : CALL section_vals_val_get(dft_section, "EPS_U_J_LOOP", r_val=dft_control%eps_u_j_loop)
485 9142 : CALL section_vals_val_get(dft_section, "MAX_MTLR_LOOP", i_val=dft_control%max_mtlr_iter)
486 : CALL section_vals_val_get(dft_section, "MTLR_REFERENCE_SCF", &
487 9142 : explicit=dft_control%mtlr_reference_scf_explicit)
488 9142 : IF (dft_control%mtlr_reference_scf_explicit) THEN
489 : CALL section_vals_val_get(dft_section, "MTLR_REFERENCE_SCF", &
490 10 : l_val=dft_control%mtlr_reference_scf)
491 : END IF
492 :
493 : ! Smearing in use
494 9142 : dft_control%smear = .FALSE.
495 :
496 : ! Surface dipole correction
497 9142 : dft_control%correct_surf_dip = .FALSE.
498 9142 : CALL section_vals_val_get(dft_section, "SURFACE_DIPOLE_CORRECTION", l_val=dft_control%correct_surf_dip)
499 9142 : CALL section_vals_val_get(dft_section, "SURF_DIP_DIR", i_val=dft_control%dir_surf_dip)
500 9142 : dft_control%pos_dir_surf_dip = -1.0_dp
501 9142 : CALL section_vals_val_get(dft_section, "SURF_DIP_POS", r_val=dft_control%pos_dir_surf_dip)
502 : ! another logical variable, surf_dip_correct_switch, is introduced for
503 : ! implementation of "SURF_DIP_SWITCH" [SGh]
504 9142 : dft_control%switch_surf_dip = .FALSE.
505 9142 : dft_control%surf_dip_correct_switch = dft_control%correct_surf_dip
506 9142 : CALL section_vals_val_get(dft_section, "SURF_DIP_SWITCH", l_val=dft_control%switch_surf_dip)
507 9142 : dft_control%correct_el_density_dip = .FALSE.
508 9142 : CALL section_vals_val_get(dft_section, "CORE_CORR_DIP", l_val=dft_control%correct_el_density_dip)
509 9142 : IF (dft_control%correct_el_density_dip) THEN
510 4 : IF (dft_control%correct_surf_dip) THEN
511 : ! Do nothing, move on
512 : ELSE
513 0 : dft_control%correct_el_density_dip = .FALSE.
514 0 : CPWARN("CORE_CORR_DIP keyword is activated only if SURFACE_DIPOLE_CORRECTION is TRUE")
515 : END IF
516 : END IF
517 :
518 : CALL section_vals_val_get(dft_section, "BASIS_SET_FILE_NAME", &
519 9142 : c_val=basis_set_file_name)
520 : CALL section_vals_val_get(dft_section, "POTENTIAL_FILE_NAME", &
521 9142 : c_val=potential_file_name)
522 :
523 : ! Read the input section
524 9142 : tmp_section => section_vals_get_subs_vals(dft_section, "sic")
525 : CALL section_vals_val_get(tmp_section, "SIC_METHOD", &
526 9142 : i_val=dft_control%sic_method_id)
527 : CALL section_vals_val_get(tmp_section, "ORBITAL_SET", &
528 9142 : i_val=dft_control%sic_list_id)
529 : CALL section_vals_val_get(tmp_section, "SIC_SCALING_A", &
530 9142 : r_val=dft_control%sic_scaling_a)
531 : CALL section_vals_val_get(tmp_section, "SIC_SCALING_B", &
532 9142 : r_val=dft_control%sic_scaling_b)
533 :
534 9142 : do_rtp = .FALSE.
535 9142 : tmp_section => section_vals_get_subs_vals(dft_section, "REAL_TIME_PROPAGATION")
536 9142 : CALL section_vals_get(tmp_section, explicit=is_present)
537 9142 : IF (is_present) THEN
538 324 : CALL read_rtp_section(dft_control, tmp_section)
539 324 : do_rtp = .TRUE.
540 : END IF
541 :
542 : ! Read the input section
543 9142 : tmp_section => section_vals_get_subs_vals(dft_section, "XAS")
544 9142 : CALL section_vals_get(tmp_section, explicit=dft_control%do_xas_calculation)
545 9142 : IF (dft_control%do_xas_calculation) THEN
546 : ! Override with section parameter
547 : CALL section_vals_val_get(tmp_section, "_SECTION_PARAMETERS_", &
548 42 : l_val=dft_control%do_xas_calculation)
549 : END IF
550 :
551 9142 : tmp_section => section_vals_get_subs_vals(dft_section, "XAS_TDP")
552 9142 : CALL section_vals_get(tmp_section, explicit=dft_control%do_xas_tdp_calculation)
553 9142 : IF (dft_control%do_xas_tdp_calculation) THEN
554 : ! Override with section parameter
555 : CALL section_vals_val_get(tmp_section, "_SECTION_PARAMETERS_", &
556 52 : l_val=dft_control%do_xas_tdp_calculation)
557 : END IF
558 :
559 : ! Read the finite field input section
560 9142 : dft_control%apply_efield = .FALSE.
561 9142 : dft_control%apply_efield_field = .FALSE. !this is for RTP
562 9142 : dft_control%apply_vector_potential = .FALSE. !this is for RTP
563 9142 : tmp_section => section_vals_get_subs_vals(dft_section, "EFIELD")
564 9142 : CALL section_vals_get(tmp_section, n_repetition=nrep, explicit=is_present)
565 9142 : IF (is_present) THEN
566 1368 : ALLOCATE (dft_control%efield_fields(nrep))
567 342 : CALL read_efield_sections(dft_control, tmp_section, cell)
568 342 : IF (do_rtp) THEN
569 30 : IF (.NOT. dft_control%rtp_control%velocity_gauge) THEN
570 20 : dft_control%apply_efield_field = .TRUE.
571 : ELSE
572 10 : dft_control%apply_vector_potential = .TRUE.
573 : ! Use this input value of vector potential to (re)start RTP
574 40 : dft_control%rtp_control%vec_pot = dft_control%efield_fields(1)%efield%vec_pot_initial
575 : END IF
576 : ELSE
577 312 : dft_control%apply_efield = .TRUE.
578 312 : CPASSERT(nrep == 1)
579 : END IF
580 : END IF
581 :
582 : ! Now, can try to guess polarisation in rtp
583 9142 : IF (do_rtp) THEN
584 : ! tmp_section => section_vals_get_subs_vals(dft_section, "REAL_TIME_PROPAGATION%PRINT%POLARIZABILITY")
585 : ! CALL section_vals_get(tmp_section, explicit=is_present)
586 : local_moment_possible = (dft_control%rtp_control%rtp_method == rtp_method_bse .OR. &
587 : dft_control%rtp_control%rtp_method == rtp_method_bse_linearized) .OR. &
588 324 : ((.NOT. dft_control%rtp_control%periodic) .AND. dft_control%rtp_control%linear_scaling)
589 90 : IF (local_moment_possible .AND. (.NOT. ASSOCIATED(dft_control%rtp_control%print_pol_elements))) THEN
590 90 : tmp_section => section_vals_get_subs_vals(dft_section, "REAL_TIME_PROPAGATION")
591 : CALL guess_pol_elements(dft_control, &
592 90 : dft_control%rtp_control%print_pol_elements)
593 : END IF
594 : END IF
595 :
596 : ! Read the finite field input section for periodic fields
597 9142 : tmp_section => section_vals_get_subs_vals(dft_section, "PERIODIC_EFIELD")
598 9142 : CALL section_vals_get(tmp_section, explicit=dft_control%apply_period_efield)
599 9142 : IF (dft_control%apply_period_efield) THEN
600 532 : ALLOCATE (dft_control%period_efield)
601 76 : CALL section_vals_val_get(tmp_section, "POLARISATION", r_vals=pol)
602 532 : dft_control%period_efield%polarisation(1:3) = pol(1:3)
603 76 : IF (PRESENT(cell)) THEN
604 76 : IF (ASSOCIATED(cell)) THEN
605 76 : CALL cell_transform_input_cartesian(cell, dft_control%period_efield%polarisation(1:3))
606 : END IF
607 : END IF
608 76 : CALL section_vals_val_get(tmp_section, "D_FILTER", r_vals=pol)
609 532 : dft_control%period_efield%d_filter(1:3) = pol(1:3)
610 76 : IF (PRESENT(cell)) THEN
611 76 : IF (ASSOCIATED(cell)) CALL cell_transform_input_cartesian(cell, dft_control%period_efield%d_filter(1:3))
612 : END IF
613 : CALL section_vals_val_get(tmp_section, "INTENSITY", &
614 76 : r_val=dft_control%period_efield%strength)
615 76 : dft_control%period_efield%displacement_field = .FALSE.
616 : CALL section_vals_val_get(tmp_section, "DISPLACEMENT_FIELD", &
617 76 : l_val=dft_control%period_efield%displacement_field)
618 :
619 76 : CALL section_vals_val_get(tmp_section, "INTENSITY_LIST", r_vals=pol)
620 :
621 76 : CALL section_vals_val_get(tmp_section, "INTENSITIES_FILE_NAME", c_val=intensities_file_name)
622 :
623 76 : IF (SIZE(pol) > 1 .OR. pol(1) /= 0.0_dp) THEN
624 : ! if INTENSITY_LIST is present, INTENSITY and INTENSITIES_FILE_NAME must not be present
625 2 : IF (dft_control%period_efield%strength /= 0.0_dp .OR. intensities_file_name /= "") THEN
626 : CALL cp_abort(__LOCATION__, "[PERIODIC FIELD] Only one of INTENSITY, INTENSITY_LIST "// &
627 0 : "or INTENSITIES_FILE_NAME can be specified.")
628 : END IF
629 :
630 6 : ALLOCATE (dft_control%period_efield%strength_list(SIZE(pol)))
631 50 : dft_control%period_efield%strength_list(1:SIZE(pol)) = pol(1:SIZE(pol))
632 : END IF
633 :
634 76 : IF (intensities_file_name /= "") THEN
635 : ! if INTENSITIES_FILE_NAME is present, INTENSITY must not be present
636 2 : IF (dft_control%period_efield%strength /= 0.0_dp) THEN
637 : CALL cp_abort(__LOCATION__, "[PERIODIC FIELD] Only one of INTENSITY, INTENSITY_LIST "// &
638 0 : "or INTENSITIES_FILE_NAME can be specified.")
639 : END IF
640 :
641 2 : CALL parser_create(parser, intensities_file_name)
642 :
643 2 : nrep = 0
644 24 : DO WHILE (.TRUE.)
645 26 : CALL parser_read_line(parser, 1, at_end)
646 26 : IF (at_end) EXIT
647 24 : nrep = nrep + 1
648 : END DO
649 :
650 2 : IF (nrep == 0) THEN
651 0 : CPABORT("[PERIODIC FIELD] No intensities found in INTENSITIES_FILE_NAME")
652 : END IF
653 :
654 6 : ALLOCATE (dft_control%period_efield%strength_list(nrep))
655 :
656 2 : CALL parser_reset(parser)
657 26 : DO irep = 1, nrep
658 24 : CALL parser_read_line(parser, 1)
659 26 : READ (parser%input_line, *) dft_control%period_efield%strength_list(irep)
660 : END DO
661 :
662 4 : CALL parser_release(parser)
663 : END IF
664 :
665 : CALL section_vals_val_get(tmp_section, "START_FRAME", &
666 76 : i_val=dft_control%period_efield%start_frame)
667 : CALL section_vals_val_get(tmp_section, "END_FRAME", &
668 76 : i_val=dft_control%period_efield%end_frame)
669 :
670 76 : IF (dft_control%period_efield%end_frame /= -1) THEN
671 : ! check if valid bounds are given
672 : ! if an end frame is given, the number of active frames must be a
673 : ! multiple of the number of intensities
674 4 : IF (dft_control%period_efield%start_frame > dft_control%period_efield%end_frame) THEN
675 0 : CPABORT("[PERIODIC FIELD] START_FRAME > END_FRAME")
676 4 : ELSE IF (dft_control%period_efield%start_frame < 1) THEN
677 0 : CPABORT("[PERIODIC FIELD] START_FRAME < 1")
678 4 : ELSE IF (MOD(dft_control%period_efield%end_frame - &
679 : dft_control%period_efield%start_frame + 1, SIZE(pol)) /= 0) THEN
680 : CALL cp_abort(__LOCATION__, &
681 0 : "[PERIODIC FIELD] Number of active frames must be a multiple of the number of intensities")
682 : END IF
683 : END IF
684 :
685 : ! periodic fields don't work with RTP
686 76 : IF (do_rtp) THEN
687 : CALL cp_abort(__LOCATION__, &
688 : "Periodic efield cannot be used with RTP. When restarting a "// &
689 : "run with periodic efield, set RESTART_RTP under &EXT_RESTART "// &
690 0 : "section to .FALSE. explicitly if RESTART_DEFAULT is .TRUE.")
691 : END IF
692 76 : IF (dft_control%period_efield%displacement_field) THEN
693 16 : CALL cite_reference(Stengel2009)
694 : ELSE
695 60 : CALL cite_reference(Souza2002)
696 60 : CALL cite_reference(Umari2002)
697 : END IF
698 : END IF
699 :
700 : ! Read the external potential input section
701 9142 : tmp_section => section_vals_get_subs_vals(dft_section, "EXTERNAL_POTENTIAL")
702 9142 : CALL section_vals_get(tmp_section, explicit=dft_control%apply_external_potential)
703 9142 : IF (dft_control%apply_external_potential) THEN
704 16 : CALL expot_control_create(dft_control%expot_control)
705 : CALL section_vals_val_get(tmp_section, "READ_FROM_CUBE", &
706 16 : l_val=dft_control%expot_control%read_from_cube)
707 : CALL section_vals_val_get(tmp_section, "STATIC", &
708 16 : l_val=dft_control%expot_control%static)
709 : CALL section_vals_val_get(tmp_section, "SCALING_FACTOR", &
710 16 : r_val=dft_control%expot_control%scaling_factor)
711 : ! External potential using Maxwell equation
712 16 : maxwell_section => section_vals_get_subs_vals(tmp_section, "MAXWELL")
713 16 : CALL section_vals_get(maxwell_section, explicit=is_present)
714 16 : IF (is_present) THEN
715 0 : dft_control%expot_control%maxwell_solver = .TRUE.
716 0 : CALL maxwell_control_create(dft_control%maxwell_control)
717 : ! read the input values from Maxwell section
718 : CALL section_vals_val_get(maxwell_section, "TEST_REAL", &
719 0 : r_val=dft_control%maxwell_control%real_test)
720 : CALL section_vals_val_get(maxwell_section, "TEST_INTEGER", &
721 0 : i_val=dft_control%maxwell_control%int_test)
722 : CALL section_vals_val_get(maxwell_section, "TEST_LOGICAL", &
723 0 : l_val=dft_control%maxwell_control%log_test)
724 : ELSE
725 16 : dft_control%expot_control%maxwell_solver = .FALSE.
726 : END IF
727 : END IF
728 :
729 : ! Read the SCCS input section if present
730 9142 : sccs_section => section_vals_get_subs_vals(dft_section, "SCCS")
731 9142 : CALL section_vals_get(sccs_section, explicit=is_present)
732 9142 : IF (is_present) THEN
733 : ! Check section parameter if SCCS is activated
734 : CALL section_vals_val_get(sccs_section, "_SECTION_PARAMETERS_", &
735 12 : l_val=dft_control%do_sccs)
736 12 : IF (dft_control%do_sccs) THEN
737 12 : ALLOCATE (dft_control%sccs_control)
738 : CALL section_vals_val_get(sccs_section, "RELATIVE_PERMITTIVITY", &
739 12 : r_val=dft_control%sccs_control%epsilon_solvent)
740 : CALL section_vals_val_get(sccs_section, "ALPHA", &
741 12 : r_val=dft_control%sccs_control%alpha_solvent)
742 : CALL section_vals_val_get(sccs_section, "BETA", &
743 12 : r_val=dft_control%sccs_control%beta_solvent)
744 : CALL section_vals_val_get(sccs_section, "DELTA_RHO", &
745 12 : r_val=dft_control%sccs_control%delta_rho)
746 : CALL section_vals_val_get(sccs_section, "DERIVATIVE_METHOD", &
747 12 : i_val=dft_control%sccs_control%derivative_method)
748 : CALL section_vals_val_get(sccs_section, "METHOD", &
749 12 : i_val=dft_control%sccs_control%method_id)
750 : CALL section_vals_val_get(sccs_section, "EPS_SCCS", &
751 12 : r_val=dft_control%sccs_control%eps_sccs)
752 : CALL section_vals_val_get(sccs_section, "EPS_SCF", &
753 12 : r_val=dft_control%sccs_control%eps_scf)
754 : CALL section_vals_val_get(sccs_section, "GAMMA", &
755 12 : r_val=dft_control%sccs_control%gamma_solvent)
756 : CALL section_vals_val_get(sccs_section, "MAX_ITER", &
757 12 : i_val=dft_control%sccs_control%max_iter)
758 : CALL section_vals_val_get(sccs_section, "MIXING", &
759 12 : r_val=dft_control%sccs_control%mixing)
760 22 : SELECT CASE (dft_control%sccs_control%method_id)
761 : CASE (sccs_andreussi)
762 10 : tmp_section => section_vals_get_subs_vals(sccs_section, "ANDREUSSI")
763 : CALL section_vals_val_get(tmp_section, "RHO_MAX", &
764 10 : r_val=dft_control%sccs_control%rho_max)
765 : CALL section_vals_val_get(tmp_section, "RHO_MIN", &
766 10 : r_val=dft_control%sccs_control%rho_min)
767 10 : IF (dft_control%sccs_control%rho_max < dft_control%sccs_control%rho_min) THEN
768 : CALL cp_abort(__LOCATION__, &
769 : "The SCCS parameter RHO_MAX is smaller than RHO_MIN. "// &
770 0 : "Please, check your input!")
771 : END IF
772 10 : CALL cite_reference(Andreussi2012)
773 : CASE (sccs_fattebert_gygi)
774 2 : tmp_section => section_vals_get_subs_vals(sccs_section, "FATTEBERT-GYGI")
775 : CALL section_vals_val_get(tmp_section, "BETA", &
776 2 : r_val=dft_control%sccs_control%beta)
777 2 : IF (dft_control%sccs_control%beta < 0.5_dp) THEN
778 : CALL cp_abort(__LOCATION__, &
779 : "A value smaller than 0.5 for the SCCS parameter beta "// &
780 0 : "causes numerical problems. Please, check your input!")
781 : END IF
782 : CALL section_vals_val_get(tmp_section, "RHO_ZERO", &
783 2 : r_val=dft_control%sccs_control%rho_zero)
784 2 : CALL cite_reference(Fattebert2002)
785 : CASE (sccs_saa_andreussi)
786 0 : tmp_section => section_vals_get_subs_vals(sccs_section, "SAA_ANDREUSSI")
787 0 : CALL section_vals_get(tmp_section, explicit=is_present)
788 0 : IF (.NOT. is_present) THEN
789 : CALL cp_abort(__LOCATION__, &
790 : "SCCS method SAA_ANDREUSSI requires the "// &
791 0 : "SAA_ANDREUSSI section.")
792 : END IF
793 0 : IF (.NOT. ALL(cell%perd == 1)) THEN
794 : CALL cp_abort(__LOCATION__, &
795 : "SCCS method SAA_ANDREUSSI is only implemented for "// &
796 0 : "3D periodic calculations.")
797 : END IF
798 : CALL section_vals_val_get(tmp_section, "RHO_MAX", &
799 0 : r_val=dft_control%sccs_control%rho_max)
800 : CALL section_vals_val_get(tmp_section, "RHO_MIN", &
801 0 : r_val=dft_control%sccs_control%rho_min)
802 0 : IF (dft_control%sccs_control%rho_max < dft_control%sccs_control%rho_min) THEN
803 : CALL cp_abort(__LOCATION__, &
804 : "The SCCS parameter RHO_MAX is smaller than RHO_MIN. "// &
805 0 : "Please, check your input!")
806 : END IF
807 : CALL section_vals_val_get(tmp_section, "F0", &
808 0 : r_val=dft_control%sccs_control%f0)
809 : CALL section_vals_val_get(tmp_section, "DELTA_ETA", &
810 0 : r_val=dft_control%sccs_control%delta_eta)
811 : CALL section_vals_val_get(tmp_section, "DELTA_ZETA", &
812 0 : r_val=dft_control%sccs_control%delta_zeta)
813 : CALL section_vals_val_get(tmp_section, "R_SOLV", &
814 0 : r_val=dft_control%sccs_control%r_solv)
815 : CALL section_vals_val_get(tmp_section, "ALPHA_ZETA", &
816 0 : r_val=dft_control%sccs_control%alpha_zeta)
817 0 : CALL cite_reference(Andreussi2019)
818 0 : CALL cite_reference(Chai2025a)
819 : CASE DEFAULT
820 12 : CPABORT("Invalid SCCS model specified. Please, check your input!")
821 : END SELECT
822 12 : CALL cite_reference(Yin2017)
823 : END IF
824 : END IF
825 :
826 : ! Read the planar counter charge section
827 9142 : tmp_section => section_vals_get_subs_vals(dft_section, "PLANAR_COUNTER_CHARGE")
828 9142 : CALL section_vals_get(tmp_section, explicit=is_present)
829 9142 : IF (is_present) THEN
830 : ! Check section parameter if planar counter charge is activated
831 : CALL section_vals_val_get(tmp_section, "_SECTION_PARAMETERS_", &
832 6 : l_val=dft_control%do_pcc)
833 6 : IF (dft_control%do_pcc) THEN
834 6 : ALLOCATE (dft_control%pcc_control)
835 : CALL section_vals_val_get(tmp_section, "DIST_EDGE", &
836 6 : r_val=dft_control%pcc_control%dist_edge)
837 : CALL section_vals_val_get(tmp_section, "GAU_C", &
838 6 : r_val=dft_control%pcc_control%gau_c)
839 : CALL section_vals_val_get(tmp_section, "PARALLEL_PLANE", &
840 6 : i_val=dft_control%pcc_control%surf_normal)
841 : END IF
842 : END IF
843 :
844 : ! Read the planar averaged Hartree potential section
845 9142 : tmp_section => section_vals_get_subs_vals(dft_section, "PLANAR_AVERAGED_V_HARTREE")
846 9142 : CALL section_vals_get(tmp_section, explicit=is_present)
847 9142 : IF (is_present) THEN
848 : ! Check section parameter if the planar-averaged potential is activated
849 : CALL section_vals_val_get(tmp_section, "_SECTION_PARAMETERS_", &
850 2 : l_val=dft_control%do_paep)
851 2 : IF (dft_control%do_paep) THEN
852 2 : ALLOCATE (dft_control%paep_control)
853 : CALL section_vals_val_get(tmp_section, "PARALLEL_PLANE", &
854 2 : i_val=dft_control%paep_control%surf_normal)
855 : END IF
856 : END IF
857 :
858 : ! ZMP added input sections
859 : ! Read the external density input section
860 9142 : tmp_section => section_vals_get_subs_vals(dft_section, "EXTERNAL_DENSITY")
861 9142 : CALL section_vals_get(tmp_section, explicit=dft_control%apply_external_density)
862 :
863 : ! Read the external vxc input section
864 9142 : tmp_section => section_vals_get_subs_vals(dft_section, "EXTERNAL_VXC")
865 9142 : CALL section_vals_get(tmp_section, explicit=dft_control%apply_external_vxc)
866 :
867 9142 : tmp_section => section_vals_get_subs_vals(dft_section, "LOCALIZE")
868 9142 : CALL section_vals_val_get(tmp_section, "EACH", i_val=dft_control%localize_each)
869 :
870 : ! SMEAGOL interface
871 9142 : tmp_section => section_vals_get_subs_vals(dft_section, "SMEAGOL")
872 9142 : CALL read_smeagol_control(dft_control%smeagol_control, tmp_section)
873 :
874 27426 : END SUBROUTINE read_dft_control
875 :
876 : ! **************************************************************************************************
877 : !> \brief Reads the input and stores in the rixs_control_type
878 : !> \param rixs_control ...
879 : !> \param rixs_section ...
880 : !> \param qs_control ...
881 : ! **************************************************************************************************
882 32 : SUBROUTINE read_rixs_control(rixs_control, rixs_section, qs_control)
883 : TYPE(rixs_control_type), POINTER :: rixs_control
884 : TYPE(section_vals_type), POINTER :: rixs_section
885 : TYPE(qs_control_type), POINTER :: qs_control
886 :
887 : TYPE(section_vals_type), POINTER :: td_section, xas_section
888 :
889 32 : CALL section_vals_val_get(rixs_section, "_SECTION_PARAMETERS_", l_val=rixs_control%enabled)
890 :
891 32 : CALL section_vals_val_get(rixs_section, "CORE_STATES", i_val=rixs_control%core_states)
892 32 : CALL section_vals_val_get(rixs_section, "VALENCE_STATES", i_val=rixs_control%valence_states)
893 :
894 32 : td_section => section_vals_get_subs_vals(rixs_section, "TDDFPT")
895 32 : CALL read_tddfpt2_control(rixs_control%tddfpt2_control, td_section, qs_control)
896 :
897 32 : xas_section => section_vals_get_subs_vals(rixs_section, "XAS_TDP")
898 32 : CALL read_xas_tdp_control(rixs_control%xas_tdp_control, xas_section)
899 :
900 32 : END SUBROUTINE read_rixs_control
901 :
902 : ! **************************************************************************************************
903 : !> \brief ...
904 : !> \param qs_control ...
905 : !> \param dft_section ...
906 : ! **************************************************************************************************
907 9142 : SUBROUTINE read_mgrid_section(qs_control, dft_section)
908 :
909 : TYPE(qs_control_type), INTENT(INOUT) :: qs_control
910 : TYPE(section_vals_type), POINTER :: dft_section
911 :
912 : CHARACTER(len=*), PARAMETER :: routineN = 'read_mgrid_section'
913 :
914 : INTEGER :: handle, igrid_level, interp_kind, &
915 : ngrid_level
916 : LOGICAL :: explicit, multigrid_set
917 : REAL(dp) :: cutoff
918 9142 : REAL(dp), DIMENSION(:), POINTER :: cutofflist
919 : TYPE(section_vals_type), POINTER :: interp_section, mgrid_section
920 :
921 9142 : CALL timeset(routineN, handle)
922 :
923 9142 : NULLIFY (interp_section, mgrid_section, cutofflist)
924 9142 : mgrid_section => section_vals_get_subs_vals(dft_section, "MGRID")
925 9142 : interp_section => section_vals_get_subs_vals(mgrid_section, "INTERPOLATOR")
926 :
927 9142 : CALL section_vals_val_get(mgrid_section, "NGRIDS", i_val=ngrid_level)
928 9142 : CALL section_vals_val_get(mgrid_section, "MULTIGRID_SET", l_val=multigrid_set)
929 9142 : CALL section_vals_val_get(mgrid_section, "CUTOFF", r_val=cutoff)
930 9142 : CALL section_vals_val_get(mgrid_section, "PROGRESSION_FACTOR", r_val=qs_control%progression_factor)
931 9142 : CALL section_vals_val_get(mgrid_section, "COMMENSURATE", l_val=qs_control%commensurate_mgrids)
932 9142 : CALL section_vals_val_get(interp_section, "KIND", i_val=interp_kind)
933 9142 : IF (interp_kind /= pw_interp) qs_control%commensurate_mgrids = .TRUE.
934 9142 : CALL section_vals_val_get(mgrid_section, "REALSPACE", l_val=qs_control%realspace_mgrids)
935 9142 : CALL section_vals_val_get(mgrid_section, "REL_CUTOFF", r_val=qs_control%relative_cutoff)
936 : CALL section_vals_val_get(mgrid_section, "SKIP_LOAD_BALANCE_DISTRIBUTED", &
937 9142 : l_val=qs_control%skip_load_balance_distributed)
938 :
939 : ! For SE and DFTB possibly override with new defaults
940 9142 : IF (qs_control%semi_empirical .OR. qs_control%dftb .OR. qs_control%xtb) THEN
941 2534 : ngrid_level = 1
942 2534 : multigrid_set = .FALSE.
943 : ! Override default cutoff value unless user specified an explicit argument..
944 2534 : CALL section_vals_val_get(mgrid_section, "CUTOFF", explicit=explicit, r_val=cutoff)
945 2534 : IF (.NOT. explicit) cutoff = 1.0_dp
946 : END IF
947 :
948 27426 : ALLOCATE (qs_control%e_cutoff(ngrid_level))
949 9142 : qs_control%cutoff = cutoff
950 :
951 9142 : IF (multigrid_set) THEN
952 : ! Read the values from input
953 4 : IF (qs_control%commensurate_mgrids) THEN
954 0 : CPABORT("Do not specify cutoffs for the commensurate grids (NYI)")
955 : END IF
956 :
957 4 : CALL section_vals_val_get(mgrid_section, "MULTIGRID_CUTOFF", r_vals=cutofflist)
958 4 : IF (ASSOCIATED(cutofflist)) THEN
959 4 : IF (SIZE(cutofflist, 1) /= ngrid_level) THEN
960 0 : CPABORT("Number of multi-grids requested and number of cutoff values do not match")
961 : END IF
962 20 : DO igrid_level = 1, ngrid_level
963 20 : qs_control%e_cutoff(igrid_level) = cutofflist(igrid_level)
964 : END DO
965 : END IF
966 : ! set cutoff to smallest value in multgrid available with >= cutoff
967 20 : DO igrid_level = ngrid_level, 1, -1
968 16 : IF (qs_control%cutoff <= qs_control%e_cutoff(igrid_level)) THEN
969 0 : qs_control%cutoff = qs_control%e_cutoff(igrid_level)
970 0 : EXIT
971 : END IF
972 : ! set largest grid value to cutoff
973 20 : IF (igrid_level == 1) THEN
974 4 : qs_control%cutoff = qs_control%e_cutoff(1)
975 : END IF
976 : END DO
977 : ELSE
978 9138 : IF (qs_control%commensurate_mgrids) qs_control%progression_factor = 4.0_dp
979 9138 : qs_control%e_cutoff(1) = qs_control%cutoff
980 28652 : DO igrid_level = 2, ngrid_level
981 : qs_control%e_cutoff(igrid_level) = qs_control%e_cutoff(igrid_level - 1)/ &
982 28652 : qs_control%progression_factor
983 : END DO
984 : END IF
985 : ! check that multigrids are ordered
986 28668 : DO igrid_level = 2, ngrid_level
987 28668 : IF (qs_control%e_cutoff(igrid_level) > qs_control%e_cutoff(igrid_level - 1)) THEN
988 0 : CPABORT("The cutoff values for the multi-grids are not ordered from large to small")
989 19526 : ELSE IF (qs_control%e_cutoff(igrid_level) == qs_control%e_cutoff(igrid_level - 1)) THEN
990 0 : CPABORT("The same cutoff value was specified for two multi-grids")
991 : END IF
992 : END DO
993 9142 : CALL timestop(handle)
994 27426 : END SUBROUTINE read_mgrid_section
995 :
996 : ! **************************************************************************************************
997 : !> \brief ...
998 : !> \param qs_control ...
999 : !> \param qs_section ...
1000 : !> \param cell optional cell used to transform Cartesian input vectors
1001 : ! **************************************************************************************************
1002 146272 : SUBROUTINE read_qs_section(qs_control, qs_section, cell)
1003 :
1004 : TYPE(qs_control_type), INTENT(INOUT) :: qs_control
1005 : TYPE(section_vals_type), POINTER :: qs_section
1006 : TYPE(cell_type), OPTIONAL, POINTER :: cell
1007 :
1008 : CHARACTER(len=*), PARAMETER :: routineN = 'read_qs_section'
1009 :
1010 : CHARACTER(LEN=2) :: element_symbol
1011 : CHARACTER(LEN=default_string_length) :: cval
1012 : CHARACTER(LEN=default_string_length), &
1013 9142 : DIMENSION(:), POINTER :: clist
1014 : INTEGER :: handle, itmp, j, jj, k, n_rep, n_var, &
1015 : ngauss, ngp, nrep, znum
1016 9142 : INTEGER, DIMENSION(:), POINTER :: tmplist
1017 : LOGICAL :: dftb_scc_mixer_explicit, dftb_tblite_mixer_explicit, explicit, &
1018 : tblite_reference_cli, tblite_reference_cli_section, tblite_section_active, was_present, &
1019 : xtb_scc_mixer_explicit, xtb_tblite_mixer_explicit
1020 : REAL(dp) :: tmp, tmpsqrt, value
1021 9142 : REAL(dp), POINTER :: scal(:)
1022 : TYPE(section_vals_type), POINTER :: cdft_control_section, ddapc_restraint_section, &
1023 : dftb_parameter, dftb_section, dftb_tblite_mixer, eeq_section, genpot_section, &
1024 : lri_optbas_section, mull_section, nonbonded_section, s2_restraint_section, se_section, &
1025 : xtb_parameter, xtb_section, xtb_tblite, xtb_tblite_mixer, xtb_tblite_ref_cli
1026 :
1027 9142 : CALL timeset(routineN, handle)
1028 :
1029 9142 : was_present = .FALSE.
1030 9142 : NULLIFY (mull_section, ddapc_restraint_section, s2_restraint_section, &
1031 9142 : se_section, dftb_section, xtb_section, dftb_parameter, xtb_parameter, lri_optbas_section, &
1032 9142 : cdft_control_section, genpot_section, eeq_section, dftb_tblite_mixer, &
1033 9142 : xtb_tblite_mixer, xtb_tblite_ref_cli)
1034 :
1035 9142 : mull_section => section_vals_get_subs_vals(qs_section, "MULLIKEN_RESTRAINT")
1036 9142 : ddapc_restraint_section => section_vals_get_subs_vals(qs_section, "DDAPC_RESTRAINT")
1037 9142 : s2_restraint_section => section_vals_get_subs_vals(qs_section, "S2_RESTRAINT")
1038 9142 : se_section => section_vals_get_subs_vals(qs_section, "SE")
1039 9142 : dftb_section => section_vals_get_subs_vals(qs_section, "DFTB")
1040 9142 : xtb_section => section_vals_get_subs_vals(qs_section, "xTB")
1041 9142 : dftb_parameter => section_vals_get_subs_vals(dftb_section, "PARAMETER")
1042 9142 : dftb_tblite_mixer => section_vals_get_subs_vals(dftb_section, "TBLITE_MIXER")
1043 9142 : xtb_parameter => section_vals_get_subs_vals(xtb_section, "PARAMETER")
1044 9142 : eeq_section => section_vals_get_subs_vals(xtb_section, "EEQ")
1045 9142 : lri_optbas_section => section_vals_get_subs_vals(qs_section, "OPTIMIZE_LRI_BASIS")
1046 9142 : cdft_control_section => section_vals_get_subs_vals(qs_section, "CDFT")
1047 9142 : nonbonded_section => section_vals_get_subs_vals(xtb_section, "NONBONDED")
1048 9142 : genpot_section => section_vals_get_subs_vals(nonbonded_section, "GENPOT")
1049 9142 : xtb_tblite_mixer => section_vals_get_subs_vals(xtb_section, "TBLITE_MIXER")
1050 9142 : xtb_tblite => section_vals_get_subs_vals(xtb_section, "TBLITE")
1051 9142 : xtb_tblite_ref_cli => section_vals_get_subs_vals(xtb_tblite, "REFERENCE_CLI")
1052 :
1053 : ! Setup all defaults values and overwrite input parameters
1054 : ! EPS_DEFAULT should set the target accuracy in the total energy (~per electron) or a closely related value
1055 9142 : CALL section_vals_val_get(qs_section, "EPS_DEFAULT", r_val=value)
1056 9142 : tmpsqrt = SQRT(value) ! a trick to work around a NAG 5.1 optimizer bug
1057 :
1058 : ! random choice ?
1059 9142 : qs_control%eps_core_charge = value/100.0_dp
1060 : ! correct if all Gaussians would have the same radius (overlap will be smaller than eps_pgf_orb**2).
1061 : ! Can be significantly in error if not... requires fully new screening/pairlist procedures
1062 9142 : qs_control%eps_pgf_orb = tmpsqrt
1063 9142 : qs_control%eps_kg_orb = qs_control%eps_pgf_orb
1064 : ! consistent since also a kind of overlap
1065 9142 : qs_control%eps_ppnl = qs_control%eps_pgf_orb/100.0_dp
1066 : ! accuracy is basically set by the overlap, this sets an empirical shift
1067 9142 : qs_control%eps_ppl = 1.0E-2_dp
1068 : !
1069 9142 : qs_control%gapw_control%eps_cpc = value
1070 : ! expexted error in the density
1071 9142 : qs_control%eps_rho_gspace = value
1072 9142 : qs_control%eps_rho_rspace = value
1073 : ! error in the gradient, can be the sqrt of the error in the energy, ignored if map_consistent
1074 9142 : qs_control%eps_gvg_rspace = tmpsqrt
1075 : !
1076 9142 : CALL section_vals_val_get(qs_section, "EPS_CORE_CHARGE", n_rep_val=n_rep)
1077 9142 : IF (n_rep /= 0) THEN
1078 0 : CALL section_vals_val_get(qs_section, "EPS_CORE_CHARGE", r_val=qs_control%eps_core_charge)
1079 : END IF
1080 9142 : CALL section_vals_val_get(qs_section, "EPS_GVG_RSPACE", n_rep_val=n_rep)
1081 9142 : IF (n_rep /= 0) THEN
1082 164 : CALL section_vals_val_get(qs_section, "EPS_GVG_RSPACE", r_val=qs_control%eps_gvg_rspace)
1083 : END IF
1084 9142 : CALL section_vals_val_get(qs_section, "EPS_PGF_ORB", n_rep_val=n_rep)
1085 9142 : IF (n_rep /= 0) THEN
1086 672 : CALL section_vals_val_get(qs_section, "EPS_PGF_ORB", r_val=qs_control%eps_pgf_orb)
1087 : END IF
1088 9142 : CALL section_vals_val_get(qs_section, "EPS_KG_ORB", n_rep_val=n_rep)
1089 9142 : IF (n_rep /= 0) THEN
1090 66 : CALL section_vals_val_get(qs_section, "EPS_KG_ORB", r_val=tmp)
1091 66 : qs_control%eps_kg_orb = SQRT(tmp)
1092 : END IF
1093 9142 : CALL section_vals_val_get(qs_section, "EPS_PPL", n_rep_val=n_rep)
1094 9142 : IF (n_rep /= 0) THEN
1095 9142 : CALL section_vals_val_get(qs_section, "EPS_PPL", r_val=qs_control%eps_ppl)
1096 : END IF
1097 9142 : CALL section_vals_val_get(qs_section, "EPS_PPNL", n_rep_val=n_rep)
1098 9142 : IF (n_rep /= 0) THEN
1099 8 : CALL section_vals_val_get(qs_section, "EPS_PPNL", r_val=qs_control%eps_ppnl)
1100 : END IF
1101 9142 : CALL section_vals_val_get(qs_section, "EPS_RHO", n_rep_val=n_rep)
1102 9142 : IF (n_rep /= 0) THEN
1103 30 : CALL section_vals_val_get(qs_section, "EPS_RHO", r_val=qs_control%eps_rho_gspace)
1104 30 : qs_control%eps_rho_rspace = qs_control%eps_rho_gspace
1105 : END IF
1106 9142 : CALL section_vals_val_get(qs_section, "EPS_RHO_RSPACE", n_rep_val=n_rep)
1107 9142 : IF (n_rep /= 0) THEN
1108 2 : CALL section_vals_val_get(qs_section, "EPS_RHO_RSPACE", r_val=qs_control%eps_rho_rspace)
1109 : END IF
1110 9142 : CALL section_vals_val_get(qs_section, "EPS_RHO_GSPACE", n_rep_val=n_rep)
1111 9142 : IF (n_rep /= 0) THEN
1112 2 : CALL section_vals_val_get(qs_section, "EPS_RHO_GSPACE", r_val=qs_control%eps_rho_gspace)
1113 : END IF
1114 9142 : CALL section_vals_val_get(qs_section, "EPS_FILTER_MATRIX", n_rep_val=n_rep)
1115 9142 : IF (n_rep /= 0) THEN
1116 9142 : CALL section_vals_val_get(qs_section, "EPS_FILTER_MATRIX", r_val=qs_control%eps_filter_matrix)
1117 : END IF
1118 9142 : CALL section_vals_val_get(qs_section, "EPS_CPC", n_rep_val=n_rep)
1119 9142 : IF (n_rep /= 0) THEN
1120 0 : CALL section_vals_val_get(qs_section, "EPS_CPC", r_val=qs_control%gapw_control%eps_cpc)
1121 : END IF
1122 :
1123 9142 : CALL section_vals_val_get(qs_section, "EPSFIT", r_val=qs_control%gapw_control%eps_fit)
1124 9142 : CALL section_vals_val_get(qs_section, "EPSISO", r_val=qs_control%gapw_control%eps_iso)
1125 9142 : CALL section_vals_val_get(qs_section, "EPSSVD", r_val=qs_control%gapw_control%eps_svd)
1126 9142 : CALL section_vals_val_get(qs_section, "EPSRHO0", r_val=qs_control%gapw_control%eps_Vrho0)
1127 9142 : CALL section_vals_val_get(qs_section, "ALPHA0_HARD", r_val=qs_control%gapw_control%alpha0_hard)
1128 9142 : qs_control%gapw_control%alpha0_hard_from_input = .FALSE.
1129 9142 : IF (qs_control%gapw_control%alpha0_hard /= 0.0_dp) qs_control%gapw_control%alpha0_hard_from_input = .TRUE.
1130 9142 : CALL section_vals_val_get(qs_section, "FORCE_PAW", l_val=qs_control%gapw_control%force_paw)
1131 9142 : CALL section_vals_val_get(qs_section, "MAX_RAD_LOCAL", r_val=qs_control%gapw_control%max_rad_local)
1132 :
1133 9142 : CALL section_vals_val_get(qs_section, "MIN_PAIR_LIST_RADIUS", r_val=qs_control%pairlist_radius)
1134 :
1135 9142 : CALL section_vals_val_get(qs_section, "LS_SCF", l_val=qs_control%do_ls_scf)
1136 9142 : CALL section_vals_val_get(qs_section, "ALMO_SCF", l_val=qs_control%do_almo_scf)
1137 9142 : CALL section_vals_val_get(qs_section, "KG_METHOD", l_val=qs_control%do_kg)
1138 :
1139 : ! Logicals
1140 9142 : CALL section_vals_val_get(qs_section, "REF_EMBED_SUBSYS", l_val=qs_control%ref_embed_subsys)
1141 9142 : CALL section_vals_val_get(qs_section, "CLUSTER_EMBED_SUBSYS", l_val=qs_control%cluster_embed_subsys)
1142 9142 : CALL section_vals_val_get(qs_section, "HIGH_LEVEL_EMBED_SUBSYS", l_val=qs_control%high_level_embed_subsys)
1143 9142 : CALL section_vals_val_get(qs_section, "DFET_EMBEDDED", l_val=qs_control%dfet_embedded)
1144 9142 : CALL section_vals_val_get(qs_section, "DMFET_EMBEDDED", l_val=qs_control%dmfet_embedded)
1145 :
1146 : ! Integers gapw
1147 9142 : CALL section_vals_val_get(qs_section, "LMAXN1", i_val=qs_control%gapw_control%lmax_sphere)
1148 9142 : CALL section_vals_val_get(qs_section, "LMAXN0", i_val=qs_control%gapw_control%lmax_rho0)
1149 9142 : CALL section_vals_val_get(qs_section, "LADDN0", i_val=qs_control%gapw_control%ladd_rho0)
1150 9142 : CALL section_vals_val_get(qs_section, "QUADRATURE", i_val=qs_control%gapw_control%quadrature)
1151 : ! GAPW 1c basis
1152 9142 : CALL section_vals_val_get(qs_section, "GAPW_1C_BASIS", i_val=qs_control%gapw_control%basis_1c)
1153 9142 : IF (qs_control%gapw_control%basis_1c /= gapw_1c_orb) THEN
1154 106 : qs_control%gapw_control%eps_svd = MAX(qs_control%gapw_control%eps_svd, 1.E-12_dp)
1155 : END IF
1156 : ! GAPW accurate integration
1157 9142 : CALL section_vals_val_get(qs_section, "GAPW_ACCURATE_XCINT", l_val=qs_control%gapw_control%accurate_xcint)
1158 9142 : CALL section_vals_val_get(qs_section, "ALPHA_WEIGHTS", r_val=qs_control%gapw_control%aweights)
1159 9142 : CALL section_vals_val_get(qs_section, "ORDER_WEIGHTS", i_val=qs_control%gapw_control%oweights)
1160 :
1161 : ! Integers grids
1162 9142 : CALL section_vals_val_get(qs_section, "PW_GRID", i_val=itmp)
1163 0 : SELECT CASE (itmp)
1164 : CASE (do_pwgrid_spherical)
1165 0 : qs_control%pw_grid_opt%spherical = .TRUE.
1166 0 : qs_control%pw_grid_opt%fullspace = .FALSE.
1167 : CASE (do_pwgrid_ns_fullspace)
1168 9142 : qs_control%pw_grid_opt%spherical = .FALSE.
1169 9142 : qs_control%pw_grid_opt%fullspace = .TRUE.
1170 : CASE (do_pwgrid_ns_halfspace)
1171 0 : qs_control%pw_grid_opt%spherical = .FALSE.
1172 9142 : qs_control%pw_grid_opt%fullspace = .FALSE.
1173 : END SELECT
1174 :
1175 : ! Method for PPL calculation
1176 9142 : CALL section_vals_val_get(qs_section, "CORE_PPL", i_val=itmp)
1177 9142 : qs_control%do_ppl_method = itmp
1178 :
1179 9142 : CALL section_vals_val_get(qs_section, "PW_GRID_LAYOUT", i_vals=tmplist)
1180 27426 : qs_control%pw_grid_opt%distribution_layout = tmplist
1181 9142 : CALL section_vals_val_get(qs_section, "PW_GRID_BLOCKED", i_val=qs_control%pw_grid_opt%blocked)
1182 :
1183 : !Integers extrapolation
1184 9142 : CALL section_vals_val_get(qs_section, "EXTRAPOLATION", i_val=qs_control%wf_interpolation_method_nr)
1185 9142 : CALL section_vals_val_get(qs_section, "EXTRAPOLATION_ORDER", i_val=qs_control%wf_extrapolation_order)
1186 :
1187 : !Method
1188 9142 : CALL section_vals_val_get(qs_section, "METHOD", i_val=qs_control%method_id)
1189 9142 : qs_control%gapw = .FALSE.
1190 9142 : qs_control%gapw_xc = .FALSE.
1191 9142 : qs_control%gpw = .FALSE.
1192 9142 : qs_control%pao = .FALSE.
1193 9142 : qs_control%dftb = .FALSE.
1194 9142 : qs_control%xtb = .FALSE.
1195 9142 : qs_control%semi_empirical = .FALSE.
1196 9142 : qs_control%ofgpw = .FALSE.
1197 9142 : qs_control%lrigpw = .FALSE.
1198 9142 : qs_control%rigpw = .FALSE.
1199 10416 : SELECT CASE (qs_control%method_id)
1200 : CASE (do_method_gapw)
1201 1274 : CALL cite_reference(Lippert1999)
1202 1274 : CALL cite_reference(Krack2000)
1203 1274 : qs_control%gapw = .TRUE.
1204 : CASE (do_method_gapw_xc)
1205 184 : qs_control%gapw_xc = .TRUE.
1206 : CASE (do_method_gpw)
1207 5106 : CALL cite_reference(Lippert1997)
1208 5106 : CALL cite_reference(VandeVondele2005a)
1209 5106 : qs_control%gpw = .TRUE.
1210 : CASE (do_method_ofgpw)
1211 0 : qs_control%ofgpw = .TRUE.
1212 : CASE (do_method_lrigpw)
1213 42 : qs_control%lrigpw = .TRUE.
1214 : CASE (do_method_rigpw)
1215 2 : qs_control%rigpw = .TRUE.
1216 : CASE (do_method_dftb)
1217 298 : qs_control%dftb = .TRUE.
1218 298 : CALL cite_reference(Porezag1995)
1219 298 : CALL cite_reference(Seifert1996)
1220 : CASE (do_method_xtb)
1221 1236 : qs_control%xtb = .TRUE.
1222 1236 : CALL cite_reference(Grimme2017)
1223 1236 : CALL cite_reference(Pracht2019)
1224 : CASE (do_method_mndo)
1225 52 : CALL cite_reference(Dewar1977)
1226 52 : qs_control%semi_empirical = .TRUE.
1227 : CASE (do_method_am1)
1228 112 : CALL cite_reference(Dewar1985)
1229 112 : qs_control%semi_empirical = .TRUE.
1230 : CASE (do_method_pm3)
1231 48 : CALL cite_reference(Stewart1989)
1232 48 : qs_control%semi_empirical = .TRUE.
1233 : CASE (do_method_pnnl)
1234 14 : CALL cite_reference(Schenter2008)
1235 14 : qs_control%semi_empirical = .TRUE.
1236 : CASE (do_method_pm6)
1237 754 : CALL cite_reference(Stewart2007)
1238 754 : qs_control%semi_empirical = .TRUE.
1239 : CASE (do_method_pm6fm)
1240 0 : CALL cite_reference(VanVoorhis2015)
1241 0 : qs_control%semi_empirical = .TRUE.
1242 : CASE (do_method_pdg)
1243 2 : CALL cite_reference(Repasky2002)
1244 2 : qs_control%semi_empirical = .TRUE.
1245 : CASE (do_method_rm1)
1246 2 : CALL cite_reference(Rocha2006)
1247 2 : qs_control%semi_empirical = .TRUE.
1248 : CASE (do_method_mndod)
1249 16 : CALL cite_reference(Dewar1977)
1250 16 : CALL cite_reference(Thiel1992)
1251 9158 : qs_control%semi_empirical = .TRUE.
1252 : END SELECT
1253 :
1254 9142 : CALL section_vals_get(mull_section, explicit=qs_control%mulliken_restraint)
1255 :
1256 9142 : IF (qs_control%mulliken_restraint) THEN
1257 2 : CALL section_vals_val_get(mull_section, "STRENGTH", r_val=qs_control%mulliken_restraint_control%strength)
1258 2 : CALL section_vals_val_get(mull_section, "TARGET", r_val=qs_control%mulliken_restraint_control%target)
1259 2 : CALL section_vals_val_get(mull_section, "ATOMS", n_rep_val=n_rep)
1260 2 : jj = 0
1261 4 : DO k = 1, n_rep
1262 2 : CALL section_vals_val_get(mull_section, "ATOMS", i_rep_val=k, i_vals=tmplist)
1263 4 : jj = jj + SIZE(tmplist)
1264 : END DO
1265 2 : qs_control%mulliken_restraint_control%natoms = jj
1266 2 : IF (qs_control%mulliken_restraint_control%natoms < 1) THEN
1267 0 : CPABORT("Need at least 1 atom to use mulliken constraints")
1268 : END IF
1269 6 : ALLOCATE (qs_control%mulliken_restraint_control%atoms(qs_control%mulliken_restraint_control%natoms))
1270 2 : jj = 0
1271 6 : DO k = 1, n_rep
1272 2 : CALL section_vals_val_get(mull_section, "ATOMS", i_rep_val=k, i_vals=tmplist)
1273 6 : DO j = 1, SIZE(tmplist)
1274 2 : jj = jj + 1
1275 4 : qs_control%mulliken_restraint_control%atoms(jj) = tmplist(j)
1276 : END DO
1277 : END DO
1278 : END IF
1279 9142 : CALL section_vals_get(ddapc_restraint_section, n_repetition=nrep, explicit=qs_control%ddapc_restraint)
1280 9142 : IF (qs_control%ddapc_restraint) THEN
1281 60 : ALLOCATE (qs_control%ddapc_restraint_control(nrep))
1282 14 : CALL read_ddapc_section(qs_control, qs_section=qs_section)
1283 14 : qs_control%ddapc_restraint_is_spin = .FALSE.
1284 14 : qs_control%ddapc_explicit_potential = .FALSE.
1285 : END IF
1286 :
1287 9142 : CALL section_vals_get(s2_restraint_section, explicit=qs_control%s2_restraint)
1288 9142 : IF (qs_control%s2_restraint) THEN
1289 : CALL section_vals_val_get(s2_restraint_section, "STRENGTH", &
1290 0 : r_val=qs_control%s2_restraint_control%strength)
1291 : CALL section_vals_val_get(s2_restraint_section, "TARGET", &
1292 0 : r_val=qs_control%s2_restraint_control%target)
1293 : CALL section_vals_val_get(s2_restraint_section, "FUNCTIONAL_FORM", &
1294 0 : i_val=qs_control%s2_restraint_control%functional_form)
1295 : END IF
1296 :
1297 9142 : CALL section_vals_get(cdft_control_section, explicit=qs_control%cdft)
1298 9142 : IF (qs_control%cdft) THEN
1299 298 : CALL read_cdft_control_section(qs_control, cdft_control_section)
1300 : END IF
1301 :
1302 : ! Semi-empirical code
1303 9142 : IF (qs_control%semi_empirical) THEN
1304 : CALL section_vals_val_get(se_section, "ORTHOGONAL_BASIS", &
1305 1000 : l_val=qs_control%se_control%orthogonal_basis)
1306 : CALL section_vals_val_get(se_section, "DELTA", &
1307 1000 : r_val=qs_control%se_control%delta)
1308 : CALL section_vals_val_get(se_section, "ANALYTICAL_GRADIENTS", &
1309 1000 : l_val=qs_control%se_control%analytical_gradients)
1310 : CALL section_vals_val_get(se_section, "FORCE_KDSO-D_EXCHANGE", &
1311 1000 : l_val=qs_control%se_control%force_kdsod_EX)
1312 : ! Integral Screening
1313 : CALL section_vals_val_get(se_section, "INTEGRAL_SCREENING", &
1314 1000 : i_val=qs_control%se_control%integral_screening)
1315 1000 : IF (qs_control%method_id == do_method_pnnl) THEN
1316 14 : IF (qs_control%se_control%integral_screening /= do_se_IS_slater) THEN
1317 : CALL cp_warn(__LOCATION__, &
1318 : "PNNL semi-empirical parameterization supports only the Slater type "// &
1319 0 : "integral scheme. Revert to Slater and continue the calculation.")
1320 : END IF
1321 14 : qs_control%se_control%integral_screening = do_se_IS_slater
1322 : END IF
1323 : ! Global Arrays variable
1324 : CALL section_vals_val_get(se_section, "GA%NCELLS", &
1325 1000 : i_val=qs_control%se_control%ga_ncells)
1326 : ! Long-Range correction
1327 : CALL section_vals_val_get(se_section, "LR_CORRECTION%CUTOFF", &
1328 1000 : r_val=qs_control%se_control%cutoff_lrc)
1329 1000 : qs_control%se_control%taper_lrc = qs_control%se_control%cutoff_lrc
1330 : CALL section_vals_val_get(se_section, "LR_CORRECTION%RC_TAPER", &
1331 1000 : explicit=explicit)
1332 1000 : IF (explicit) THEN
1333 : CALL section_vals_val_get(se_section, "LR_CORRECTION%RC_TAPER", &
1334 0 : r_val=qs_control%se_control%taper_lrc)
1335 : END IF
1336 : CALL section_vals_val_get(se_section, "LR_CORRECTION%RC_RANGE", &
1337 1000 : r_val=qs_control%se_control%range_lrc)
1338 : ! Coulomb
1339 : CALL section_vals_val_get(se_section, "COULOMB%CUTOFF", &
1340 1000 : r_val=qs_control%se_control%cutoff_cou)
1341 1000 : qs_control%se_control%taper_cou = qs_control%se_control%cutoff_cou
1342 : CALL section_vals_val_get(se_section, "COULOMB%RC_TAPER", &
1343 1000 : explicit=explicit)
1344 1000 : IF (explicit) THEN
1345 : CALL section_vals_val_get(se_section, "COULOMB%RC_TAPER", &
1346 0 : r_val=qs_control%se_control%taper_cou)
1347 : END IF
1348 : CALL section_vals_val_get(se_section, "COULOMB%RC_RANGE", &
1349 1000 : r_val=qs_control%se_control%range_cou)
1350 : ! Exchange
1351 : CALL section_vals_val_get(se_section, "EXCHANGE%CUTOFF", &
1352 1000 : r_val=qs_control%se_control%cutoff_exc)
1353 1000 : qs_control%se_control%taper_exc = qs_control%se_control%cutoff_exc
1354 : CALL section_vals_val_get(se_section, "EXCHANGE%RC_TAPER", &
1355 1000 : explicit=explicit)
1356 1000 : IF (explicit) THEN
1357 : CALL section_vals_val_get(se_section, "EXCHANGE%RC_TAPER", &
1358 38 : r_val=qs_control%se_control%taper_exc)
1359 : END IF
1360 : CALL section_vals_val_get(se_section, "EXCHANGE%RC_RANGE", &
1361 1000 : r_val=qs_control%se_control%range_exc)
1362 : ! Screening (only if the integral scheme is of dumped type)
1363 1000 : IF (qs_control%se_control%integral_screening == do_se_IS_kdso_d) THEN
1364 : CALL section_vals_val_get(se_section, "SCREENING%RC_TAPER", &
1365 14 : r_val=qs_control%se_control%taper_scr)
1366 : CALL section_vals_val_get(se_section, "SCREENING%RC_RANGE", &
1367 14 : r_val=qs_control%se_control%range_scr)
1368 : END IF
1369 : ! Periodic Type Calculation
1370 : CALL section_vals_val_get(se_section, "PERIODIC", &
1371 1000 : i_val=qs_control%se_control%periodic_type)
1372 1968 : SELECT CASE (qs_control%se_control%periodic_type)
1373 : CASE (do_se_lr_none)
1374 968 : qs_control%se_control%do_ewald = .FALSE.
1375 968 : qs_control%se_control%do_ewald_r3 = .FALSE.
1376 968 : qs_control%se_control%do_ewald_gks = .FALSE.
1377 : CASE (do_se_lr_ewald)
1378 30 : qs_control%se_control%do_ewald = .TRUE.
1379 30 : qs_control%se_control%do_ewald_r3 = .FALSE.
1380 30 : qs_control%se_control%do_ewald_gks = .FALSE.
1381 : CASE (do_se_lr_ewald_gks)
1382 2 : qs_control%se_control%do_ewald = .FALSE.
1383 2 : qs_control%se_control%do_ewald_r3 = .FALSE.
1384 2 : qs_control%se_control%do_ewald_gks = .TRUE.
1385 2 : IF (qs_control%method_id /= do_method_pnnl) THEN
1386 : CALL cp_abort(__LOCATION__, &
1387 : "A periodic semi-empirical calculation was requested with a long-range "// &
1388 : "summation on the single integral evaluation. This scheme is supported "// &
1389 0 : "only by the PNNL parameterization.")
1390 : END IF
1391 : CASE (do_se_lr_ewald_r3)
1392 0 : qs_control%se_control%do_ewald = .TRUE.
1393 0 : qs_control%se_control%do_ewald_r3 = .TRUE.
1394 0 : qs_control%se_control%do_ewald_gks = .FALSE.
1395 1000 : IF (qs_control%se_control%integral_screening /= do_se_IS_kdso) THEN
1396 : CALL cp_abort(__LOCATION__, &
1397 : "A periodic semi-empirical calculation was requested with a long-range "// &
1398 : "summation for the slowly convergent part 1/R^3, which is not congruent "// &
1399 : "with the integral screening chosen. The only integral screening supported "// &
1400 0 : "by this periodic type calculation is the standard Klopman-Dewar-Sabelli-Ohno.")
1401 : END IF
1402 : END SELECT
1403 :
1404 : ! dispersion pair potentials
1405 : CALL section_vals_val_get(se_section, "DISPERSION", &
1406 1000 : l_val=qs_control%se_control%dispersion)
1407 : CALL section_vals_val_get(se_section, "DISPERSION_RADIUS", &
1408 1000 : r_val=qs_control%se_control%rcdisp)
1409 : CALL section_vals_val_get(se_section, "COORDINATION_CUTOFF", &
1410 1000 : r_val=qs_control%se_control%epscn)
1411 1000 : CALL section_vals_val_get(se_section, "D3_SCALING", r_vals=scal)
1412 1000 : qs_control%se_control%sd3(1) = scal(1)
1413 1000 : qs_control%se_control%sd3(2) = scal(2)
1414 1000 : qs_control%se_control%sd3(3) = scal(3)
1415 : CALL section_vals_val_get(se_section, "DISPERSION_PARAMETER_FILE", &
1416 1000 : c_val=qs_control%se_control%dispersion_parameter_file)
1417 :
1418 : ! Stop the execution for non-implemented features
1419 1000 : IF (qs_control%se_control%periodic_type == do_se_lr_ewald_r3) THEN
1420 0 : CPABORT("EWALD_R3 not implemented yet!")
1421 : END IF
1422 :
1423 : IF (qs_control%method_id == do_method_mndo .OR. &
1424 : qs_control%method_id == do_method_am1 .OR. &
1425 : qs_control%method_id == do_method_mndod .OR. &
1426 : qs_control%method_id == do_method_pdg .OR. &
1427 : qs_control%method_id == do_method_pm3 .OR. &
1428 : qs_control%method_id == do_method_pm6 .OR. &
1429 : qs_control%method_id == do_method_pm6fm .OR. &
1430 1000 : qs_control%method_id == do_method_pnnl .OR. &
1431 : qs_control%method_id == do_method_rm1) THEN
1432 1000 : qs_control%se_control%orthogonal_basis = .TRUE.
1433 : END IF
1434 : END IF
1435 :
1436 : ! DFTB code
1437 9142 : IF (qs_control%dftb) THEN
1438 : CALL section_vals_val_get(dftb_section, "ORTHOGONAL_BASIS", &
1439 298 : l_val=qs_control%dftb_control%orthogonal_basis)
1440 : CALL section_vals_val_get(dftb_section, "SELF_CONSISTENT", &
1441 298 : l_val=qs_control%dftb_control%self_consistent)
1442 : CALL section_vals_val_get(dftb_section, "DISPERSION", &
1443 298 : l_val=qs_control%dftb_control%dispersion)
1444 : CALL section_vals_val_get(dftb_section, "DIAGONAL_DFTB3", &
1445 298 : l_val=qs_control%dftb_control%dftb3_diagonal)
1446 : CALL section_vals_val_get(dftb_section, "HB_SR_GAMMA", &
1447 298 : l_val=qs_control%dftb_control%hb_sr_damp)
1448 : CALL section_vals_val_get(dftb_section, "SCC_MIXER", &
1449 298 : explicit=dftb_scc_mixer_explicit)
1450 : CALL section_vals_val_get(dftb_section, "SCC_MIXER", &
1451 298 : i_val=qs_control%dftb_control%tblite_scc_mixer)
1452 298 : CALL section_vals_get(dftb_tblite_mixer, explicit=dftb_tblite_mixer_explicit)
1453 : CALL read_tblite_mixer_section(dftb_tblite_mixer, &
1454 : qs_control%dftb_control%tblite_mixer_iterations, &
1455 : qs_control%dftb_control%tblite_mixer_memory, &
1456 : qs_control%dftb_control%tblite_mixer_solver, &
1457 : qs_control%dftb_control%tblite_mixer_damping, &
1458 : qs_control%dftb_control%tblite_mixer_omega0, &
1459 : qs_control%dftb_control%tblite_mixer_min_weight, &
1460 : qs_control%dftb_control%tblite_mixer_max_weight, &
1461 : qs_control%dftb_control%tblite_mixer_weight_factor, &
1462 298 : "DFTB/TBLITE_MIXER")
1463 298 : IF (qs_control%do_ls_scf) THEN
1464 44 : IF (dftb_scc_mixer_explicit .AND. &
1465 : qs_control%dftb_control%tblite_scc_mixer /= tblite_scc_mixer_none) THEN
1466 : CALL cp_warn(__LOCATION__, &
1467 : "DFTB/SCC_MIXER is reset to NONE with QS/LS_SCF; LS_SCF optimizes "// &
1468 2 : "the density matrix directly.")
1469 : END IF
1470 44 : IF (dftb_tblite_mixer_explicit) THEN
1471 : CALL cp_warn(__LOCATION__, &
1472 : "DFTB/TBLITE_MIXER settings are ignored with QS/LS_SCF; LS_SCF controls "// &
1473 0 : "the density-matrix optimization.")
1474 : END IF
1475 44 : qs_control%dftb_control%tblite_scc_mixer = tblite_scc_mixer_none
1476 : END IF
1477 298 : IF (qs_control%dftb_control%tblite_mixer_damping <= 0.0_dp) THEN
1478 0 : CPABORT("DFTB/TBLITE_MIXER/DAMPING must be positive")
1479 : END IF
1480 : CALL section_vals_val_get(dftb_section, "EPS_DISP", &
1481 298 : r_val=qs_control%dftb_control%eps_disp)
1482 298 : CALL section_vals_val_get(dftb_section, "DO_EWALD", explicit=explicit)
1483 298 : IF (explicit) THEN
1484 : CALL section_vals_val_get(dftb_section, "DO_EWALD", &
1485 206 : l_val=qs_control%dftb_control%do_ewald)
1486 : ELSE
1487 92 : qs_control%dftb_control%do_ewald = (qs_control%periodicity /= 0)
1488 : END IF
1489 : CALL section_vals_val_get(dftb_parameter, "PARAM_FILE_PATH", &
1490 298 : c_val=qs_control%dftb_control%sk_file_path)
1491 : CALL section_vals_val_get(dftb_parameter, "PARAM_FILE_NAME", &
1492 298 : c_val=qs_control%dftb_control%sk_file_list)
1493 : CALL section_vals_val_get(dftb_parameter, "HB_SR_PARAM", &
1494 298 : r_val=qs_control%dftb_control%hb_sr_para)
1495 298 : CALL section_vals_val_get(dftb_parameter, "SK_FILE", n_rep_val=n_var)
1496 644 : ALLOCATE (qs_control%dftb_control%sk_pair_list(3, n_var))
1497 394 : DO k = 1, n_var
1498 : CALL section_vals_val_get(dftb_parameter, "SK_FILE", i_rep_val=k, &
1499 96 : c_vals=clist)
1500 682 : qs_control%dftb_control%sk_pair_list(1:3, k) = clist(1:3)
1501 : END DO
1502 : ! Dispersion type
1503 : CALL section_vals_val_get(dftb_parameter, "DISPERSION_TYPE", &
1504 298 : i_val=qs_control%dftb_control%dispersion_type)
1505 : CALL section_vals_val_get(dftb_parameter, "UFF_FORCE_FIELD", &
1506 298 : c_val=qs_control%dftb_control%uff_force_field)
1507 : ! D3 Dispersion
1508 : CALL section_vals_val_get(dftb_parameter, "DISPERSION_RADIUS", &
1509 298 : r_val=qs_control%dftb_control%rcdisp)
1510 : CALL section_vals_val_get(dftb_parameter, "COORDINATION_CUTOFF", &
1511 298 : r_val=qs_control%dftb_control%epscn)
1512 : CALL section_vals_val_get(dftb_parameter, "D2_EXP_PRE", &
1513 298 : r_val=qs_control%dftb_control%exp_pre)
1514 : CALL section_vals_val_get(dftb_parameter, "D2_SCALING", &
1515 298 : r_val=qs_control%dftb_control%scaling)
1516 298 : CALL section_vals_val_get(dftb_parameter, "D3_SCALING", r_vals=scal)
1517 298 : qs_control%dftb_control%sd3(1) = scal(1)
1518 298 : qs_control%dftb_control%sd3(2) = scal(2)
1519 298 : qs_control%dftb_control%sd3(3) = scal(3)
1520 298 : CALL section_vals_val_get(dftb_parameter, "D3BJ_SCALING", r_vals=scal)
1521 298 : qs_control%dftb_control%sd3bj(1) = scal(1)
1522 298 : qs_control%dftb_control%sd3bj(2) = scal(2)
1523 298 : qs_control%dftb_control%sd3bj(3) = scal(3)
1524 298 : qs_control%dftb_control%sd3bj(4) = scal(4)
1525 : CALL section_vals_val_get(dftb_parameter, "DISPERSION_PARAMETER_FILE", &
1526 298 : c_val=qs_control%dftb_control%dispersion_parameter_file)
1527 :
1528 298 : IF (qs_control%dftb_control%dispersion) CALL cite_reference(Zhechkov2005)
1529 298 : IF (qs_control%dftb_control%self_consistent) CALL cite_reference(Elstner1998)
1530 1490 : IF (qs_control%dftb_control%hb_sr_damp) CALL cite_reference(Hu2007)
1531 : END IF
1532 :
1533 : ! xTB code
1534 9142 : IF (qs_control%xtb) THEN
1535 1236 : CALL section_vals_val_get(xtb_section, "GFN_TYPE", i_val=qs_control%xtb_control%gfn_type)
1536 1236 : CALL section_vals_val_get(xtb_tblite, "_SECTION_PARAMETERS_", l_val=tblite_section_active)
1537 1236 : qs_control%xtb_control%do_tblite = (qs_control%xtb_control%gfn_type == gfn_tblite)
1538 1236 : IF (qs_control%xtb_control%do_tblite) THEN
1539 196 : IF (.NOT. tblite_section_active) THEN
1540 0 : CPABORT("XTB/GFN_TYPE TBLITE requires an XTB/TBLITE section")
1541 : END IF
1542 : ! The CP2K-internal GFN1 defaults are still used to initialize shared xTB fields.
1543 196 : qs_control%xtb_control%gfn_type = gfn1xtb
1544 1040 : ELSE IF (tblite_section_active) THEN
1545 0 : CPABORT("The XTB/TBLITE section requires XTB/GFN_TYPE TBLITE")
1546 : END IF
1547 : CALL section_vals_val_get(xtb_section, "SCC_MIXER", &
1548 1236 : explicit=xtb_scc_mixer_explicit)
1549 : CALL section_vals_val_get(xtb_section, "SCC_MIXER", &
1550 1236 : i_val=qs_control%xtb_control%tblite_scc_mixer)
1551 1236 : CALL section_vals_get(xtb_tblite_mixer, explicit=xtb_tblite_mixer_explicit)
1552 : CALL read_tblite_mixer_section(xtb_tblite_mixer, &
1553 : qs_control%xtb_control%tblite_mixer_iterations, &
1554 : qs_control%xtb_control%tblite_mixer_memory, &
1555 : qs_control%xtb_control%tblite_mixer_solver, &
1556 : qs_control%xtb_control%tblite_mixer_damping, &
1557 : qs_control%xtb_control%tblite_mixer_omega0, &
1558 : qs_control%xtb_control%tblite_mixer_min_weight, &
1559 : qs_control%xtb_control%tblite_mixer_max_weight, &
1560 : qs_control%xtb_control%tblite_mixer_weight_factor, &
1561 1236 : "XTB/TBLITE_MIXER")
1562 1236 : IF (xtb_tblite_mixer_explicit) THEN
1563 : CALL section_vals_val_get(xtb_tblite_mixer, "DAMPING", &
1564 2 : explicit=qs_control%xtb_control%tblite_mixer_damping_explicit)
1565 : END IF
1566 1236 : IF ((.NOT. qs_control%xtb_control%do_tblite) .AND. &
1567 : qs_control%xtb_control%gfn_type == 0) THEN
1568 : IF (xtb_scc_mixer_explicit .AND. &
1569 694 : qs_control%xtb_control%tblite_scc_mixer /= tblite_scc_mixer_auto .AND. &
1570 : qs_control%xtb_control%tblite_scc_mixer /= tblite_scc_mixer_none) THEN
1571 : CALL cp_warn(__LOCATION__, &
1572 : "XTB/SCC_MIXER is reset to NONE for CP2K-internal GFN0-xTB; "// &
1573 0 : "GFN0-xTB has no SCC variables to mix.")
1574 : END IF
1575 694 : IF (xtb_tblite_mixer_explicit) THEN
1576 : CALL cp_warn(__LOCATION__, &
1577 : "XTB/TBLITE_MIXER settings are ignored for CP2K-internal GFN0-xTB; "// &
1578 0 : "GFN0-xTB has no SCC variables to mix.")
1579 : END IF
1580 694 : qs_control%xtb_control%tblite_scc_mixer = tblite_scc_mixer_none
1581 : END IF
1582 1236 : IF (qs_control%do_ls_scf) THEN
1583 36 : IF (xtb_scc_mixer_explicit .AND. &
1584 : qs_control%xtb_control%tblite_scc_mixer /= tblite_scc_mixer_none) THEN
1585 : CALL cp_warn(__LOCATION__, &
1586 : "XTB/SCC_MIXER is reset to NONE with QS/LS_SCF; LS_SCF optimizes "// &
1587 4 : "the density matrix directly.")
1588 : END IF
1589 36 : IF (xtb_tblite_mixer_explicit) THEN
1590 : CALL cp_warn(__LOCATION__, &
1591 : "XTB/TBLITE_MIXER settings are ignored with QS/LS_SCF; LS_SCF controls "// &
1592 0 : "the density-matrix optimization.")
1593 : END IF
1594 36 : qs_control%xtb_control%tblite_scc_mixer = tblite_scc_mixer_none
1595 : END IF
1596 1236 : IF (qs_control%xtb_control%tblite_mixer_damping <= 0.0_dp) THEN
1597 0 : CPABORT("XTB/TBLITE_MIXER/DAMPING must be positive")
1598 : END IF
1599 1236 : CALL section_vals_val_get(xtb_section, "DO_EWALD", explicit=explicit)
1600 1236 : IF (explicit) THEN
1601 : CALL section_vals_val_get(xtb_section, "DO_EWALD", &
1602 776 : l_val=qs_control%xtb_control%do_ewald)
1603 : ELSE
1604 460 : qs_control%xtb_control%do_ewald = (qs_control%periodicity /= 0)
1605 : END IF
1606 : ! Spin Polarisation
1607 : CALL section_vals_val_get(xtb_section, "SPIN_POLARISATION", &
1608 1236 : l_val=qs_control%xtb_control%do_spinpol)
1609 : ! vdW
1610 1236 : CALL section_vals_val_get(xtb_section, "VDW_POTENTIAL", explicit=explicit)
1611 1236 : IF (explicit) THEN
1612 684 : CALL section_vals_val_get(xtb_section, "VDW_POTENTIAL", c_val=cval)
1613 684 : CALL uppercase(cval)
1614 2 : SELECT CASE (cval)
1615 : CASE ("NONE")
1616 2 : qs_control%xtb_control%vdw_type = xtb_vdw_type_none
1617 : CASE ("DFTD3")
1618 56 : qs_control%xtb_control%vdw_type = xtb_vdw_type_d3
1619 : CASE ("DFTD4")
1620 626 : qs_control%xtb_control%vdw_type = xtb_vdw_type_d4
1621 : CASE DEFAULT
1622 684 : CPABORT("vdW type")
1623 : END SELECT
1624 : ELSE
1625 568 : SELECT CASE (qs_control%xtb_control%gfn_type)
1626 : CASE (0)
1627 16 : qs_control%xtb_control%vdw_type = xtb_vdw_type_d4
1628 : CASE (1)
1629 536 : qs_control%xtb_control%vdw_type = xtb_vdw_type_d3
1630 : CASE (2)
1631 0 : qs_control%xtb_control%vdw_type = xtb_vdw_type_d4
1632 0 : CPABORT("gfn2-xtb tbd")
1633 : CASE DEFAULT
1634 552 : CPABORT("GFN type")
1635 : END SELECT
1636 : END IF
1637 : !
1638 1236 : CALL section_vals_val_get(xtb_section, "STO_NG", i_val=ngauss)
1639 1236 : qs_control%xtb_control%sto_ng = ngauss
1640 1236 : CALL section_vals_val_get(xtb_section, "HYDROGEN_STO_NG", i_val=ngauss)
1641 1236 : qs_control%xtb_control%h_sto_ng = ngauss
1642 1236 : CALL section_vals_val_get(xtb_section, "STO_FLEX", explicit=explicit)
1643 1236 : IF (explicit) THEN
1644 : CALL section_vals_val_get(xtb_section, "STO_FLEX", &
1645 6 : l_val=qs_control%xtb_control%sto_flex)
1646 : ELSE
1647 1230 : qs_control%xtb_control%sto_flex = .FALSE.
1648 : END IF
1649 : CALL section_vals_val_get(xtb_parameter, "PARAM_FILE_PATH", &
1650 1236 : c_val=qs_control%xtb_control%parameter_file_path)
1651 1236 : CALL section_vals_val_get(xtb_parameter, "PARAM_FILE_NAME", explicit=explicit)
1652 1236 : IF (explicit) THEN
1653 : CALL section_vals_val_get(xtb_parameter, "PARAM_FILE_NAME", &
1654 0 : c_val=qs_control%xtb_control%parameter_file_name)
1655 : ELSE
1656 1930 : SELECT CASE (qs_control%xtb_control%gfn_type)
1657 : CASE (0)
1658 694 : qs_control%xtb_control%parameter_file_name = "xTB0_parameters"
1659 : CASE (1)
1660 542 : qs_control%xtb_control%parameter_file_name = "xTB1_parameters"
1661 : CASE (2)
1662 0 : CPABORT("gfn2-xtb tbd")
1663 : CASE DEFAULT
1664 1236 : CPABORT("GFN type")
1665 : END SELECT
1666 : END IF
1667 : !
1668 : CALL section_vals_val_get(xtb_parameter, "SPINPOL_PARAM_FILE_NAME", &
1669 1236 : c_val=qs_control%xtb_control%spinpol_param_file_name)
1670 : ! D3 Dispersion
1671 : CALL section_vals_val_get(xtb_parameter, "DISPERSION_RADIUS", &
1672 1236 : r_val=qs_control%xtb_control%rcdisp)
1673 : CALL section_vals_val_get(xtb_parameter, "COORDINATION_CUTOFF", &
1674 1236 : r_val=qs_control%xtb_control%epscn)
1675 1236 : CALL section_vals_val_get(xtb_parameter, "D3BJ_SCALING", explicit=explicit)
1676 1236 : IF (explicit) THEN
1677 0 : CALL section_vals_val_get(xtb_parameter, "D3BJ_SCALING", r_vals=scal)
1678 0 : qs_control%xtb_control%s6 = scal(1)
1679 0 : qs_control%xtb_control%s8 = scal(2)
1680 : ELSE
1681 1930 : SELECT CASE (qs_control%xtb_control%gfn_type)
1682 : CASE (0)
1683 694 : qs_control%xtb_control%s6 = 1.00_dp
1684 694 : qs_control%xtb_control%s8 = 2.85_dp
1685 : CASE (1)
1686 542 : qs_control%xtb_control%s6 = 1.00_dp
1687 542 : qs_control%xtb_control%s8 = 2.40_dp
1688 : CASE (2)
1689 0 : CPABORT("gfn2-xtb tbd")
1690 : CASE DEFAULT
1691 1236 : CPABORT("GFN type")
1692 : END SELECT
1693 : END IF
1694 1236 : CALL section_vals_val_get(xtb_parameter, "D3BJ_PARAM", explicit=explicit)
1695 1236 : IF (explicit) THEN
1696 0 : CALL section_vals_val_get(xtb_parameter, "D3BJ_PARAM", r_vals=scal)
1697 0 : qs_control%xtb_control%a1 = scal(1)
1698 0 : qs_control%xtb_control%a2 = scal(2)
1699 : ELSE
1700 1930 : SELECT CASE (qs_control%xtb_control%gfn_type)
1701 : CASE (0)
1702 694 : qs_control%xtb_control%a1 = 0.80_dp
1703 694 : qs_control%xtb_control%a2 = 4.60_dp
1704 : CASE (1)
1705 542 : qs_control%xtb_control%a1 = 0.63_dp
1706 542 : qs_control%xtb_control%a2 = 5.00_dp
1707 : CASE (2)
1708 0 : CPABORT("gfn2-xtb tbd")
1709 : CASE DEFAULT
1710 1236 : CPABORT("GFN type")
1711 : END SELECT
1712 : END IF
1713 : CALL section_vals_val_get(xtb_parameter, "DISPERSION_PARAMETER_FILE", &
1714 1236 : c_val=qs_control%xtb_control%dispersion_parameter_file)
1715 : ! global parameters
1716 1236 : CALL section_vals_val_get(xtb_parameter, "HUCKEL_CONSTANTS", explicit=explicit)
1717 1236 : IF (explicit) THEN
1718 0 : CALL section_vals_val_get(xtb_parameter, "HUCKEL_CONSTANTS", r_vals=scal)
1719 0 : qs_control%xtb_control%ks = scal(1)
1720 0 : qs_control%xtb_control%kp = scal(2)
1721 0 : qs_control%xtb_control%kd = scal(3)
1722 0 : qs_control%xtb_control%ksp = scal(4)
1723 0 : qs_control%xtb_control%k2sh = scal(5)
1724 0 : IF (qs_control%xtb_control%gfn_type == 0) THEN
1725 : ! enforce ksp for gfn0
1726 0 : qs_control%xtb_control%ksp = 0.5_dp*(scal(1) + scal(2))
1727 : END IF
1728 : ELSE
1729 1930 : SELECT CASE (qs_control%xtb_control%gfn_type)
1730 : CASE (0)
1731 694 : qs_control%xtb_control%ks = 2.00_dp
1732 694 : qs_control%xtb_control%kp = 2.4868_dp
1733 694 : qs_control%xtb_control%kd = 2.27_dp
1734 694 : qs_control%xtb_control%ksp = 2.2434_dp
1735 694 : qs_control%xtb_control%k2sh = 1.1241_dp
1736 : CASE (1)
1737 542 : qs_control%xtb_control%ks = 1.85_dp
1738 542 : qs_control%xtb_control%kp = 2.25_dp
1739 542 : qs_control%xtb_control%kd = 2.00_dp
1740 542 : qs_control%xtb_control%ksp = 2.08_dp
1741 542 : qs_control%xtb_control%k2sh = 2.85_dp
1742 : CASE (2)
1743 0 : CPABORT("gfn2-xtb tbd")
1744 : CASE DEFAULT
1745 1236 : CPABORT("GFN type")
1746 : END SELECT
1747 : END IF
1748 1236 : CALL section_vals_val_get(xtb_parameter, "COULOMB_CONSTANTS", explicit=explicit)
1749 1236 : IF (explicit) THEN
1750 0 : CALL section_vals_val_get(xtb_parameter, "COULOMB_CONSTANTS", r_vals=scal)
1751 0 : qs_control%xtb_control%kg = scal(1)
1752 0 : qs_control%xtb_control%kf = scal(2)
1753 : ELSE
1754 1930 : SELECT CASE (qs_control%xtb_control%gfn_type)
1755 : CASE (0)
1756 694 : qs_control%xtb_control%kg = 2.00_dp
1757 694 : qs_control%xtb_control%kf = 1.50_dp
1758 : CASE (1)
1759 542 : qs_control%xtb_control%kg = 2.00_dp
1760 542 : qs_control%xtb_control%kf = 1.50_dp
1761 : CASE (2)
1762 0 : CPABORT("gfn2-xtb tbd")
1763 : CASE DEFAULT
1764 1236 : CPABORT("GFN type")
1765 : END SELECT
1766 : END IF
1767 1236 : CALL section_vals_val_get(xtb_parameter, "CN_CONSTANTS", r_vals=scal)
1768 1236 : qs_control%xtb_control%kcns = scal(1)
1769 1236 : qs_control%xtb_control%kcnp = scal(2)
1770 1236 : qs_control%xtb_control%kcnd = scal(3)
1771 : !
1772 1236 : CALL section_vals_val_get(xtb_parameter, "EN_CONSTANTS", explicit=explicit)
1773 1236 : IF (explicit) THEN
1774 0 : CALL section_vals_val_get(xtb_parameter, "EN_CONSTANTS", r_vals=scal)
1775 0 : SELECT CASE (qs_control%xtb_control%gfn_type)
1776 : CASE (0)
1777 0 : qs_control%xtb_control%ksen = scal(1)
1778 0 : qs_control%xtb_control%kpen = scal(2)
1779 0 : qs_control%xtb_control%kden = scal(3)
1780 : CASE (1)
1781 0 : qs_control%xtb_control%ken = scal(1)
1782 : CASE (2)
1783 0 : CPABORT("gfn2-xtb tbd")
1784 : CASE DEFAULT
1785 0 : CPABORT("GFN type")
1786 : END SELECT
1787 : ELSE
1788 1930 : SELECT CASE (qs_control%xtb_control%gfn_type)
1789 : CASE (0)
1790 694 : qs_control%xtb_control%ksen = 0.006_dp
1791 694 : qs_control%xtb_control%kpen = -0.001_dp
1792 694 : qs_control%xtb_control%kden = -0.002_dp
1793 : CASE (1)
1794 542 : qs_control%xtb_control%ken = -0.007_dp
1795 : CASE (2)
1796 0 : CPABORT("gfn2-xtb tbd")
1797 : CASE DEFAULT
1798 1236 : CPABORT("GFN type")
1799 : END SELECT
1800 : END IF
1801 : ! ben
1802 1236 : CALL section_vals_val_get(xtb_parameter, "BEN_CONSTANT", r_vals=scal)
1803 1236 : qs_control%xtb_control%ben = scal(1)
1804 : ! enscale (hidden parameter in repulsion
1805 1236 : CALL section_vals_val_get(xtb_parameter, "ENSCALE", explicit=explicit)
1806 1236 : IF (explicit) THEN
1807 : CALL section_vals_val_get(xtb_parameter, "ENSCALE", &
1808 0 : r_val=qs_control%xtb_control%enscale)
1809 : ELSE
1810 1930 : SELECT CASE (qs_control%xtb_control%gfn_type)
1811 : CASE (0)
1812 694 : qs_control%xtb_control%enscale = -0.09_dp
1813 : CASE (1)
1814 542 : qs_control%xtb_control%enscale = 0._dp
1815 : CASE (2)
1816 0 : CPABORT("gfn2-xtb tbd")
1817 : CASE DEFAULT
1818 1236 : CPABORT("GFN type")
1819 : END SELECT
1820 : END IF
1821 : ! XB
1822 : CALL section_vals_val_get(xtb_section, "USE_HALOGEN_CORRECTION", &
1823 1236 : l_val=qs_control%xtb_control%xb_interaction)
1824 1236 : CALL section_vals_val_get(xtb_parameter, "HALOGEN_BINDING", r_vals=scal)
1825 1236 : qs_control%xtb_control%kxr = scal(1)
1826 1236 : qs_control%xtb_control%kx2 = scal(2)
1827 : ! NONBONDED interactions
1828 : CALL section_vals_val_get(xtb_section, "DO_NONBONDED", &
1829 1236 : l_val=qs_control%xtb_control%do_nonbonded)
1830 1236 : CALL section_vals_get(nonbonded_section, explicit=explicit)
1831 1236 : IF (explicit .AND. qs_control%xtb_control%do_nonbonded) THEN
1832 6 : CALL section_vals_get(genpot_section, explicit=explicit, n_repetition=ngp)
1833 6 : IF (explicit) THEN
1834 6 : CALL pair_potential_reallocate(qs_control%xtb_control%nonbonded, 1, ngp, gp=.TRUE.)
1835 6 : CALL read_gp_section(qs_control%xtb_control%nonbonded, genpot_section, 0)
1836 : END IF
1837 : END IF !nonbonded
1838 : CALL section_vals_val_get(xtb_section, "EPS_PAIRPOTENTIAL", &
1839 1236 : r_val=qs_control%xtb_control%eps_pair)
1840 : ! SR Coulomb
1841 1236 : CALL section_vals_val_get(xtb_parameter, "COULOMB_SR_CUT", r_vals=scal)
1842 1236 : qs_control%xtb_control%coulomb_sr_cut = scal(1)
1843 1236 : CALL section_vals_val_get(xtb_parameter, "COULOMB_SR_EPS", r_vals=scal)
1844 1236 : qs_control%xtb_control%coulomb_sr_eps = scal(1)
1845 : ! XB_radius
1846 1236 : CALL section_vals_val_get(xtb_parameter, "XB_RADIUS", r_val=qs_control%xtb_control%xb_radius)
1847 : ! Kab
1848 1236 : CALL section_vals_val_get(xtb_parameter, "KAB_PARAM", n_rep_val=n_rep)
1849 : ! Coulomb
1850 1930 : SELECT CASE (qs_control%xtb_control%gfn_type)
1851 : CASE (0)
1852 694 : qs_control%xtb_control%coulomb_interaction = .FALSE.
1853 694 : qs_control%xtb_control%coulomb_lr = .FALSE.
1854 694 : qs_control%xtb_control%tb3_interaction = .FALSE.
1855 694 : qs_control%xtb_control%check_atomic_charges = .FALSE.
1856 : CALL section_vals_val_get(xtb_section, "VARIATIONAL_DIPOLE", &
1857 694 : l_val=qs_control%xtb_control%var_dipole)
1858 : CASE (1)
1859 : ! For debugging purposes
1860 : CALL section_vals_val_get(xtb_section, "COULOMB_INTERACTION", &
1861 542 : l_val=qs_control%xtb_control%coulomb_interaction)
1862 : CALL section_vals_val_get(xtb_section, "COULOMB_LR", &
1863 542 : l_val=qs_control%xtb_control%coulomb_lr)
1864 : CALL section_vals_val_get(xtb_section, "TB3_INTERACTION", &
1865 542 : l_val=qs_control%xtb_control%tb3_interaction)
1866 : ! Check for bad atomic charges
1867 : CALL section_vals_val_get(xtb_section, "CHECK_ATOMIC_CHARGES", &
1868 542 : l_val=qs_control%xtb_control%check_atomic_charges)
1869 542 : qs_control%xtb_control%var_dipole = .FALSE.
1870 : CASE (2)
1871 0 : CPABORT("gfn2-xtb tbd")
1872 : CASE DEFAULT
1873 1236 : CPABORT("GFN type")
1874 : END SELECT
1875 1236 : qs_control%xtb_control%kab_nval = n_rep
1876 1236 : IF (n_rep > 0) THEN
1877 6 : ALLOCATE (qs_control%xtb_control%kab_param(3, n_rep))
1878 6 : ALLOCATE (qs_control%xtb_control%kab_types(2, n_rep))
1879 6 : ALLOCATE (qs_control%xtb_control%kab_vals(n_rep))
1880 4 : DO j = 1, n_rep
1881 2 : CALL section_vals_val_get(xtb_parameter, "KAB_PARAM", i_rep_val=j, c_vals=clist)
1882 2 : qs_control%xtb_control%kab_param(1, j) = clist(1)
1883 : CALL get_ptable_info(clist(1) (1:2), &
1884 2 : ielement=qs_control%xtb_control%kab_types(1, j))
1885 2 : qs_control%xtb_control%kab_param(2, j) = clist(2)
1886 : CALL get_ptable_info(clist(2) (1:2), &
1887 2 : ielement=qs_control%xtb_control%kab_types(2, j))
1888 2 : qs_control%xtb_control%kab_param(3, j) = clist(3)
1889 4 : READ (clist(3), '(F10.0)') qs_control%xtb_control%kab_vals(j)
1890 : END DO
1891 : END IF
1892 :
1893 : ! Spin Polarisation
1894 1236 : CALL section_vals_val_get(xtb_parameter, "SPIN_POL_PARAM", n_rep_val=n_rep)
1895 1236 : IF (n_rep > 0) THEN
1896 6 : ALLOCATE (qs_control%xtb_control%spinpol_type(n_rep))
1897 6 : ALLOCATE (qs_control%xtb_control%spinpol_vals(6, n_rep))
1898 8 : DO j = 1, n_rep
1899 6 : CALL section_vals_val_get(xtb_parameter, "SPIN_POL_PARAM", i_rep_val=j, c_vals=clist)
1900 6 : READ (clist(1), '(A)') cval
1901 6 : element_symbol = ADJUSTL(TRIM(cval))
1902 6 : CALL get_ptable_info(element_symbol, znum)
1903 6 : qs_control%xtb_control%spinpol_type(j) = znum
1904 6 : READ (clist(2), '(F20.8)') qs_control%xtb_control%spinpol_vals(1, j)
1905 6 : READ (clist(3), '(F20.8)') qs_control%xtb_control%spinpol_vals(2, j)
1906 6 : READ (clist(4), '(F20.8)') qs_control%xtb_control%spinpol_vals(3, j)
1907 6 : READ (clist(5), '(F20.8)') qs_control%xtb_control%spinpol_vals(4, j)
1908 6 : READ (clist(6), '(F20.8)') qs_control%xtb_control%spinpol_vals(5, j)
1909 14 : READ (clist(7), '(F20.8)') qs_control%xtb_control%spinpol_vals(6, j)
1910 : END DO
1911 : END IF
1912 :
1913 1236 : IF (qs_control%xtb_control%gfn_type == 0) THEN
1914 694 : CALL section_vals_val_get(xtb_parameter, "SRB_PARAMETER", r_vals=scal)
1915 694 : qs_control%xtb_control%ksrb = scal(1)
1916 694 : qs_control%xtb_control%esrb = scal(2)
1917 694 : qs_control%xtb_control%gscal = scal(3)
1918 694 : qs_control%xtb_control%c1srb = scal(4)
1919 694 : qs_control%xtb_control%c2srb = scal(5)
1920 694 : qs_control%xtb_control%shift = scal(6)
1921 : END IF
1922 :
1923 1236 : CALL section_vals_val_get(xtb_section, "EN_SHIFT_TYPE", c_val=cval)
1924 1236 : CALL uppercase(cval)
1925 1236 : SELECT CASE (TRIM(cval))
1926 : CASE ("SELECT")
1927 0 : qs_control%xtb_control%enshift_type = 0
1928 : CASE ("MOLECULE")
1929 1236 : qs_control%xtb_control%enshift_type = 1
1930 : CASE ("CRYSTAL")
1931 0 : qs_control%xtb_control%enshift_type = 2
1932 : CASE DEFAULT
1933 1236 : CPABORT("Unknown value for EN_SHIFT_TYPE")
1934 : END SELECT
1935 :
1936 : ! EEQ solver params
1937 1236 : CALL read_eeq_param(eeq_section, qs_control%xtb_control%eeq_sparam)
1938 : END IF
1939 :
1940 : ! Optimize LRI basis set
1941 9142 : CALL section_vals_get(lri_optbas_section, explicit=qs_control%lri_optbas)
1942 :
1943 : ! Use tblite if selected through XTB/GFN_TYPE TBLITE.
1944 9142 : IF (qs_control%xtb_control%do_tblite) THEN
1945 : CALL section_vals_val_get(xtb_tblite, "METHOD", &
1946 196 : i_val=qs_control%xtb_control%tblite_method)
1947 : CALL section_vals_val_get(xtb_tblite, "PARAM", &
1948 196 : c_val=qs_control%xtb_control%tblite_param_file)
1949 : CALL section_vals_val_get(xtb_tblite, "ACCURACY", &
1950 196 : r_val=qs_control%xtb_control%tblite_accuracy)
1951 196 : IF (qs_control%xtb_control%tblite_accuracy <= 0.0_dp) THEN
1952 0 : CPABORT("XTB/TBLITE/ACCURACY must be positive")
1953 : END IF
1954 196 : IF (qs_control%xtb_control%tblite_mixer_damping <= 0.0_dp) THEN
1955 0 : CPABORT("XTB/TBLITE_MIXER/DAMPING must be positive")
1956 : END IF
1957 196 : CALL section_vals_val_get(xtb_tblite, "REFERENCE_CLI", l_val=tblite_reference_cli)
1958 196 : CALL section_vals_get(xtb_tblite_ref_cli, explicit=tblite_reference_cli_section)
1959 196 : IF (tblite_reference_cli .AND. (.NOT. tblite_reference_cli_section)) THEN
1960 0 : CPABORT("XTB/TBLITE/REFERENCE_CLI keyword requires an XTB/TBLITE/REFERENCE_CLI section")
1961 : END IF
1962 196 : IF (tblite_reference_cli .OR. tblite_reference_cli_section) THEN
1963 2 : CALL read_xtb_reference_cli_section(xtb_tblite_ref_cli, qs_control%xtb_control%reference_cli, cell)
1964 2 : qs_control%xtb_control%reference_cli%enabled = .TRUE.
1965 : END IF
1966 196 : CALL cite_reference(Katbashev2025)
1967 : ! tblite handles periodic long-range terms internally from the CP2K cell periodicity.
1968 : ! Keep xtb_control%do_ewald as read above from XTB/DO_EWALD or SUBSYS/CELL/PERIODIC,
1969 : ! matching the DFTB and CP2K-internal xTB setup.
1970 : END IF
1971 :
1972 9142 : CALL timestop(handle)
1973 9142 : END SUBROUTINE read_qs_section
1974 :
1975 : ! **************************************************************************************************
1976 : !> \brief Read a TBLITE_MIXER section.
1977 : !> \param mixer_section input section
1978 : !> \param iterations tblite SCC iteration limit
1979 : !> \param memory Broyden history length
1980 : !> \param solver native tblite electronic solver id
1981 : !> \param damping mixer damping parameter
1982 : !> \param omega0 Broyden regularization weight
1983 : !> \param min_weight minimum dynamic Broyden weight
1984 : !> \param max_weight maximum dynamic Broyden weight
1985 : !> \param weight_factor residual-to-weight scaling factor
1986 : !> \param section_name diagnostic section name
1987 : ! **************************************************************************************************
1988 1566 : SUBROUTINE read_tblite_mixer_section(mixer_section, iterations, memory, solver, damping, omega0, min_weight, &
1989 : max_weight, weight_factor, section_name)
1990 : TYPE(section_vals_type), POINTER :: mixer_section
1991 : INTEGER, INTENT(INOUT) :: iterations, memory, solver
1992 : REAL(KIND=dp), INTENT(INOUT) :: damping, omega0, min_weight, max_weight, &
1993 : weight_factor
1994 : CHARACTER(LEN=*), INTENT(IN) :: section_name
1995 :
1996 : LOGICAL :: explicit, memory_explicit
1997 :
1998 1534 : CALL section_vals_get(mixer_section, explicit=explicit)
1999 1534 : IF (.NOT. explicit) RETURN
2000 :
2001 4 : CALL section_vals_val_get(mixer_section, "ITERATIONS", i_val=iterations)
2002 4 : CALL section_vals_val_get(mixer_section, "MEMORY", explicit=memory_explicit)
2003 4 : IF (memory_explicit) CALL section_vals_val_get(mixer_section, "MEMORY", i_val=memory)
2004 4 : IF (.NOT. memory_explicit .OR. memory == tblite_mixer_memory_inherit) THEN
2005 2 : memory = iterations
2006 : END IF
2007 4 : CALL section_vals_val_get(mixer_section, "SOLVER", i_val=solver)
2008 4 : CALL section_vals_val_get(mixer_section, "DAMPING", r_val=damping)
2009 4 : CALL section_vals_val_get(mixer_section, "OMEGA0", r_val=omega0)
2010 4 : CALL section_vals_val_get(mixer_section, "MIN_WEIGHT", r_val=min_weight)
2011 4 : CALL section_vals_val_get(mixer_section, "MAX_WEIGHT", r_val=max_weight)
2012 4 : CALL section_vals_val_get(mixer_section, "WEIGHT_FACTOR", r_val=weight_factor)
2013 :
2014 4 : IF (iterations < 1) CPABORT(TRIM(section_name)//"/ITERATIONS must be positive")
2015 4 : IF (memory < 1) CPABORT(TRIM(section_name)//"/MEMORY must be positive")
2016 4 : SELECT CASE (solver)
2017 : CASE (tblite_solver_gvd, tblite_solver_gvr)
2018 : CASE DEFAULT
2019 4 : CPABORT(TRIM(section_name)//"/SOLVER must be GVD or GVR")
2020 : END SELECT
2021 4 : IF (damping <= 0.0_dp) CPABORT(TRIM(section_name)//"/DAMPING must be positive")
2022 4 : IF (omega0 <= 0.0_dp) CPABORT(TRIM(section_name)//"/OMEGA0 must be positive")
2023 4 : IF (min_weight <= 0.0_dp) CPABORT(TRIM(section_name)//"/MIN_WEIGHT must be positive")
2024 4 : IF (max_weight <= 0.0_dp) CPABORT(TRIM(section_name)//"/MAX_WEIGHT must be positive")
2025 4 : IF (max_weight < min_weight) THEN
2026 0 : CPABORT(TRIM(section_name)//"/MAX_WEIGHT must not be smaller than MIN_WEIGHT")
2027 : END IF
2028 4 : IF (weight_factor <= 0.0_dp) CPABORT(TRIM(section_name)//"/WEIGHT_FACTOR must be positive")
2029 :
2030 : END SUBROUTINE read_tblite_mixer_section
2031 :
2032 : ! **************************************************************************************************
2033 : !> \brief Read native tblite CLI reference options.
2034 : !> \param ref_cli_section input section
2035 : !> \param ref_cli reference CLI control data
2036 : !> \param cell optional cell used to transform Cartesian input vectors
2037 : ! **************************************************************************************************
2038 2 : SUBROUTINE read_xtb_reference_cli_section(ref_cli_section, ref_cli, cell)
2039 : TYPE(section_vals_type), POINTER :: ref_cli_section
2040 : TYPE(xtb_reference_cli_type), INTENT(INOUT) :: ref_cli
2041 : TYPE(cell_type), OPTIONAL, POINTER :: cell
2042 :
2043 2 : REAL(KIND=dp), DIMENSION(:), POINTER :: efield
2044 : TYPE(section_vals_type), POINTER :: fit_section, guess_section, &
2045 : param_section, solvation_section, &
2046 : tagdiff_section
2047 :
2048 2 : CALL section_vals_val_get(ref_cli_section, "_SECTION_PARAMETERS_", l_val=ref_cli%enabled)
2049 2 : CALL section_vals_val_get(ref_cli_section, "PROGRAM_NAME", c_val=ref_cli%program_name)
2050 2 : CALL section_vals_val_get(ref_cli_section, "GUESS", i_val=ref_cli%guess)
2051 2 : CALL section_vals_val_get(ref_cli_section, "WORK_DIRECTORY", c_val=ref_cli%work_directory)
2052 2 : CALL section_vals_val_get(ref_cli_section, "PREFIX", c_val=ref_cli%prefix)
2053 2 : CALL section_vals_val_get(ref_cli_section, "INPUT_FORMAT", c_val=ref_cli%input_format)
2054 2 : CALL section_vals_val_get(ref_cli_section, "RESTART", c_val=ref_cli%restart_file)
2055 2 : CALL section_vals_val_get(ref_cli_section, "GRAD", c_val=ref_cli%grad_file)
2056 2 : CALL section_vals_val_get(ref_cli_section, "JSON", c_val=ref_cli%json_file)
2057 2 : CALL section_vals_val_get(ref_cli_section, "POST_PROCESSING", c_val=ref_cli%post_processing)
2058 : CALL section_vals_val_get(ref_cli_section, "POST_PROCESSING_OUTPUT", &
2059 2 : c_val=ref_cli%post_processing_output_file)
2060 2 : CALL section_vals_val_get(ref_cli_section, "EFIELD", explicit=ref_cli%efield_active)
2061 2 : IF (ref_cli%efield_active) THEN
2062 2 : NULLIFY (efield)
2063 2 : CALL section_vals_val_get(ref_cli_section, "EFIELD", r_vals=efield)
2064 8 : ref_cli%efield = efield(1:3)
2065 2 : IF (PRESENT(cell)) THEN
2066 2 : IF (ASSOCIATED(cell)) CALL cell_transform_input_cartesian(cell, ref_cli%efield)
2067 : END IF
2068 : END IF
2069 2 : solvation_section => section_vals_get_subs_vals(ref_cli_section, "IMPLICIT_SOLVATION")
2070 2 : CALL section_vals_get(solvation_section, explicit=ref_cli%solvation_active)
2071 2 : IF (ref_cli%solvation_active) THEN
2072 2 : CALL section_vals_val_get(solvation_section, "MODEL", i_val=ref_cli%solvation_model)
2073 2 : CALL section_vals_val_get(solvation_section, "SOLVENT", c_val=ref_cli%solvation_solvent)
2074 2 : CALL section_vals_val_get(solvation_section, "BORN_KERNEL", i_val=ref_cli%solvation_born_kernel)
2075 2 : CALL section_vals_val_get(solvation_section, "SOLUTION_STATE", i_val=ref_cli%solvation_state)
2076 2 : IF (LEN_TRIM(ref_cli%solvation_solvent) == 0) THEN
2077 0 : CPABORT("REFERENCE_CLI implicit solvation needs SOLVENT")
2078 : END IF
2079 2 : IF (ref_cli%solvation_model == tblite_cli_solvation_cpcm .AND. &
2080 : ref_cli%solvation_born_kernel /= tblite_cli_born_kernel_auto) THEN
2081 0 : CPABORT("BORN_KERNEL is invalid with MODEL CPCM")
2082 : END IF
2083 2 : IF (ref_cli%solvation_state /= tblite_cli_solution_state_gsolv) THEN
2084 2 : SELECT CASE (ref_cli%solvation_model)
2085 : CASE (tblite_cli_solvation_alpb, tblite_cli_solvation_gbsa)
2086 : ! Native tblite supports solution-state shifts for parametrized named-solvent ALPB/GBSA.
2087 : CASE (tblite_cli_solvation_gbe, tblite_cli_solvation_gb, tblite_cli_solvation_cpcm)
2088 2 : CPABORT("SOLUTION_STATE is valid only for ALPB/GBSA")
2089 : END SELECT
2090 : END IF
2091 : END IF
2092 : CALL section_vals_val_get(ref_cli_section, "ELECTRONIC_TEMPERATURE_GUESS", &
2093 2 : r_val=ref_cli%electronic_temperature_guess)
2094 2 : IF (ref_cli%electronic_temperature_guess < 0.0_dp) THEN
2095 0 : CPABORT("XTB/TBLITE/REFERENCE_CLI/ELECTRONIC_TEMPERATURE_GUESS must not be negative")
2096 : END IF
2097 2 : IF (ref_cli%electronic_temperature_guess > 0.0_dp .AND. ref_cli%guess /= tblite_guess_ceh) THEN
2098 0 : CPABORT("XTB/TBLITE/REFERENCE_CLI/ELECTRONIC_TEMPERATURE_GUESS requires GUESS CEH")
2099 : END IF
2100 2 : guess_section => section_vals_get_subs_vals(ref_cli_section, "GUESS_CLI")
2101 2 : CALL section_vals_get(guess_section, explicit=ref_cli%guess_cli%enabled)
2102 2 : IF (ref_cli%guess_cli%enabled) THEN
2103 0 : CALL section_vals_val_get(guess_section, "METHOD", i_val=ref_cli%guess_cli%method)
2104 : CALL section_vals_val_get(guess_section, "ELECTRONIC_TEMPERATURE_GUESS", &
2105 0 : r_val=ref_cli%guess_cli%electronic_temperature_guess)
2106 0 : IF (ref_cli%guess_cli%electronic_temperature_guess < 0.0_dp) THEN
2107 0 : CPABORT("REFERENCE_CLI/GUESS_CLI/ELECTRONIC_TEMPERATURE_GUESS must not be negative")
2108 : END IF
2109 0 : CALL section_vals_val_get(guess_section, "SOLVER", i_val=ref_cli%guess_cli%solver)
2110 0 : CALL section_vals_val_get(guess_section, "EFIELD", explicit=ref_cli%guess_cli%efield_active)
2111 0 : IF (ref_cli%guess_cli%efield_active) THEN
2112 0 : NULLIFY (efield)
2113 0 : CALL section_vals_val_get(guess_section, "EFIELD", r_vals=efield)
2114 0 : ref_cli%guess_cli%efield = efield(1:3)
2115 0 : IF (PRESENT(cell)) THEN
2116 0 : IF (ASSOCIATED(cell)) CALL cell_transform_input_cartesian(cell, ref_cli%guess_cli%efield)
2117 : END IF
2118 : END IF
2119 0 : CALL section_vals_val_get(guess_section, "GRAD", l_val=ref_cli%guess_cli%grad)
2120 0 : CALL section_vals_val_get(guess_section, "JSON", c_val=ref_cli%guess_cli%json_file)
2121 0 : CALL section_vals_val_get(guess_section, "INPUT_FORMAT", c_val=ref_cli%guess_cli%input_format)
2122 0 : CALL section_vals_val_get(guess_section, "INPUT_FILE", c_val=ref_cli%guess_cli%input_file)
2123 : END IF
2124 2 : param_section => section_vals_get_subs_vals(ref_cli_section, "PARAM_CLI")
2125 2 : CALL section_vals_get(param_section, explicit=ref_cli%param_cli%enabled)
2126 2 : IF (ref_cli%param_cli%enabled) THEN
2127 : CALL section_vals_val_get(param_section, "METHOD", explicit=ref_cli%param_cli%method_explicit, &
2128 0 : i_val=ref_cli%param_cli%method)
2129 0 : CALL section_vals_val_get(param_section, "OUTPUT", c_val=ref_cli%param_cli%output_file)
2130 0 : CALL section_vals_val_get(param_section, "INPUT_FILE", c_val=ref_cli%param_cli%input_file)
2131 : END IF
2132 2 : fit_section => section_vals_get_subs_vals(ref_cli_section, "FIT_CLI")
2133 2 : CALL section_vals_get(fit_section, explicit=ref_cli%fit_cli%enabled)
2134 2 : IF (ref_cli%fit_cli%enabled) THEN
2135 0 : CALL section_vals_val_get(fit_section, "PARAM_FILE", c_val=ref_cli%fit_cli%param_file)
2136 0 : CALL section_vals_val_get(fit_section, "INPUT_FILE", c_val=ref_cli%fit_cli%input_file)
2137 0 : CALL section_vals_val_get(fit_section, "DRY_RUN", l_val=ref_cli%fit_cli%dry_run)
2138 0 : CALL section_vals_val_get(fit_section, "COPY", c_val=ref_cli%fit_cli%copy_file)
2139 0 : IF (LEN_TRIM(ref_cli%fit_cli%param_file) == 0) THEN
2140 0 : CPABORT("XTB/TBLITE/REFERENCE_CLI/FIT_CLI needs PARAM_FILE")
2141 : END IF
2142 0 : IF (LEN_TRIM(ref_cli%fit_cli%input_file) == 0) THEN
2143 0 : CPABORT("XTB/TBLITE/REFERENCE_CLI/FIT_CLI needs INPUT_FILE")
2144 : END IF
2145 : END IF
2146 2 : tagdiff_section => section_vals_get_subs_vals(ref_cli_section, "TAGDIFF_CLI")
2147 2 : CALL section_vals_get(tagdiff_section, explicit=ref_cli%tagdiff_cli%enabled)
2148 2 : IF (ref_cli%tagdiff_cli%enabled) THEN
2149 0 : CALL section_vals_val_get(tagdiff_section, "ACTUAL", c_val=ref_cli%tagdiff_cli%actual_file)
2150 0 : CALL section_vals_val_get(tagdiff_section, "REFERENCE", c_val=ref_cli%tagdiff_cli%reference_file)
2151 0 : CALL section_vals_val_get(tagdiff_section, "FIT", l_val=ref_cli%tagdiff_cli%fit)
2152 0 : IF (LEN_TRIM(ref_cli%tagdiff_cli%actual_file) == 0) THEN
2153 0 : CPABORT("XTB/TBLITE/REFERENCE_CLI/TAGDIFF_CLI needs ACTUAL")
2154 : END IF
2155 0 : IF (LEN_TRIM(ref_cli%tagdiff_cli%reference_file) == 0) THEN
2156 0 : CPABORT("XTB/TBLITE/REFERENCE_CLI/TAGDIFF_CLI needs REFERENCE")
2157 : END IF
2158 : END IF
2159 2 : CALL section_vals_val_get(ref_cli_section, "KEEP_FILES", l_val=ref_cli%keep_files)
2160 2 : CALL section_vals_val_get(ref_cli_section, "ERROR_LIMIT", r_val=ref_cli%error_limit)
2161 2 : CALL section_vals_val_get(ref_cli_section, "STOP_ON_ERROR", l_val=ref_cli%stop_on_error)
2162 2 : CALL section_vals_val_get(ref_cli_section, "CHECK_ENERGY", l_val=ref_cli%check_energy)
2163 2 : CALL section_vals_val_get(ref_cli_section, "CHECK_FORCES", l_val=ref_cli%check_forces)
2164 2 : CALL section_vals_val_get(ref_cli_section, "CHECK_VIRIAL", l_val=ref_cli%check_virial)
2165 :
2166 2 : END SUBROUTINE read_xtb_reference_cli_section
2167 :
2168 : ! **************************************************************************************************
2169 : !> \brief Read TDDFPT-related input parameters.
2170 : !> \param t_control TDDFPT control parameters
2171 : !> \param t_section TDDFPT input section
2172 : !> \param qs_control Quickstep control parameters
2173 : ! **************************************************************************************************
2174 9174 : SUBROUTINE read_tddfpt2_control(t_control, t_section, qs_control)
2175 : TYPE(tddfpt2_control_type), POINTER :: t_control
2176 : TYPE(section_vals_type), POINTER :: t_section
2177 : TYPE(qs_control_type), POINTER :: qs_control
2178 :
2179 : CHARACTER(LEN=*), PARAMETER :: routineN = 'read_tddfpt2_control'
2180 :
2181 : CHARACTER(LEN=default_string_length), &
2182 9174 : DIMENSION(:), POINTER :: tmpstringlist
2183 : INTEGER :: handle, irep, isize, nrep
2184 9174 : INTEGER, ALLOCATABLE, DIMENSION(:) :: inds
2185 : LOGICAL :: do_ewald, do_exchange, expl, explicit, &
2186 : multigrid_set
2187 : REAL(KIND=dp) :: filter, fval, hfx
2188 : TYPE(section_vals_type), POINTER :: dipole_section, mgrid_section, &
2189 : soc_section, stda_section, xc_func, &
2190 : xc_section
2191 :
2192 9174 : CALL timeset(routineN, handle)
2193 :
2194 9174 : CALL section_vals_val_get(t_section, "_SECTION_PARAMETERS_", l_val=t_control%enabled)
2195 :
2196 9174 : CALL section_vals_val_get(t_section, "NSTATES", i_val=t_control%nstates)
2197 9174 : CALL section_vals_val_get(t_section, "MAX_ITER", i_val=t_control%niters)
2198 9174 : CALL section_vals_val_get(t_section, "MAX_KV", i_val=t_control%nkvs)
2199 9174 : CALL section_vals_val_get(t_section, "NLUMO", i_val=t_control%nlumo)
2200 9174 : CALL section_vals_val_get(t_section, "NPROC_STATE", i_val=t_control%nprocs)
2201 9174 : CALL section_vals_val_get(t_section, "KERNEL", i_val=t_control%kernel)
2202 9174 : CALL section_vals_val_get(t_section, "SPINFLIP", i_val=t_control%spinflip)
2203 9174 : CALL section_vals_val_get(t_section, "OE_CORR", i_val=t_control%oe_corr)
2204 9174 : CALL section_vals_val_get(t_section, "EV_SHIFT", r_val=t_control%ev_shift)
2205 9174 : CALL section_vals_val_get(t_section, "EOS_SHIFT", r_val=t_control%eos_shift)
2206 :
2207 9174 : CALL section_vals_val_get(t_section, "CONVERGENCE", r_val=t_control%conv)
2208 9174 : CALL section_vals_val_get(t_section, "MIN_AMPLITUDE", r_val=t_control%min_excitation_amplitude)
2209 9174 : CALL section_vals_val_get(t_section, "ORTHOGONAL_EPS", r_val=t_control%orthogonal_eps)
2210 :
2211 9174 : CALL section_vals_val_get(t_section, "RESTART", l_val=t_control%is_restart)
2212 9174 : CALL section_vals_val_get(t_section, "RKS_TRIPLETS", l_val=t_control%rks_triplets)
2213 9174 : CALL section_vals_val_get(t_section, "DO_LRIGPW", l_val=t_control%do_lrigpw)
2214 9174 : CALL section_vals_val_get(t_section, "DO_SMEARING", l_val=t_control%do_smearing)
2215 9174 : CALL section_vals_val_get(t_section, "DO_BSE", l_val=t_control%do_bse)
2216 9174 : CALL section_vals_val_get(t_section, "DO_BSE_W_ONLY", l_val=t_control%do_bse_w_only)
2217 9174 : CALL section_vals_val_get(t_section, "DO_BSE_GW_ONLY", l_val=t_control%do_bse_gw_only)
2218 9174 : CALL section_vals_val_get(t_section, "ADMM_KERNEL_CORRECTION_SYMMETRIC", l_val=t_control%admm_symm)
2219 9174 : CALL section_vals_val_get(t_section, "ADMM_KERNEL_XC_CORRECTION", l_val=t_control%admm_xc_correction)
2220 9174 : CALL section_vals_val_get(t_section, "EXCITON_DESCRIPTORS", l_val=t_control%do_exciton_descriptors)
2221 9174 : CALL section_vals_val_get(t_section, "DIRECTIONAL_EXCITON_DESCRIPTORS", l_val=t_control%do_directional_exciton_descriptors)
2222 : CALL section_vals_val_get(t_section, "DIRECTIONAL_EXCITON_CROSSCORRELATION", &
2223 9174 : l_val=t_control%do_directional_exciton_crosscorrelation, explicit=explicit)
2224 9174 : IF (explicit .AND. t_control%do_directional_exciton_crosscorrelation .AND. &
2225 : .NOT. t_control%do_directional_exciton_descriptors) THEN
2226 : CALL cp_warn(__LOCATION__, &
2227 0 : "DIRECTIONAL_EXCITON_CROSSCORRELATION has no effect without DIRECTIONAL_EXCITON_DESCRIPTORS.")
2228 : END IF
2229 :
2230 : ! read automatically generated auxiliary basis for LRI
2231 9174 : CALL section_vals_val_get(t_section, "AUTO_BASIS", n_rep_val=nrep)
2232 18348 : DO irep = 1, nrep
2233 9174 : CALL section_vals_val_get(t_section, "AUTO_BASIS", i_rep_val=irep, c_vals=tmpstringlist)
2234 18348 : IF (SIZE(tmpstringlist) == 2) THEN
2235 9174 : CALL uppercase(tmpstringlist(2))
2236 18348 : SELECT CASE (tmpstringlist(2))
2237 : CASE ("X")
2238 9174 : SELECT CASE (tmpstringlist(1))
2239 : CASE ("X")
2240 : ! Do nothing
2241 : CASE DEFAULT
2242 : CALL cp_abort(__LOCATION__, &
2243 : "AUTO_BASIS: the size <X> is invalid for the "// &
2244 : "type <"//TRIM(ADJUSTL(tmpstringlist(1)))//">; "// &
2245 : "use one of SMALL, MEDIUM, LARGE, HUGE for "// &
2246 : "the size. The syntax AUTO_BASIS X X is a "// &
2247 : "reserved case for using NO automatically "// &
2248 9174 : "generated basis sets.")
2249 : END SELECT
2250 : CASE ("SMALL")
2251 0 : isize = 0
2252 : CASE ("MEDIUM")
2253 0 : isize = 1
2254 : CASE ("LARGE")
2255 0 : isize = 2
2256 : CASE ("HUGE")
2257 0 : isize = 3
2258 : CASE DEFAULT
2259 9174 : CPABORT("Unknown basis size in AUTO_BASIS keyword:"//TRIM(tmpstringlist(1)))
2260 : END SELECT
2261 : !
2262 9174 : SELECT CASE (tmpstringlist(1))
2263 : CASE ("X")
2264 : CASE ("P_LRI_AUX")
2265 0 : t_control%auto_basis_p_lri_aux = isize
2266 : CASE DEFAULT
2267 9174 : CPABORT("Unknown basis type in AUTO_BASIS keyword:"//TRIM(tmpstringlist(1)))
2268 : END SELECT
2269 : ELSE
2270 : CALL cp_abort(__LOCATION__, &
2271 0 : "AUTO_BASIS keyword in &PROPERTIES &TDDFT section has a wrong number of arguments.")
2272 : END IF
2273 : END DO
2274 :
2275 9174 : IF (t_control%conv < 0) THEN
2276 0 : t_control%conv = ABS(t_control%conv)
2277 : END IF
2278 :
2279 : ! DIPOLE_MOMENTS subsection
2280 9174 : dipole_section => section_vals_get_subs_vals(t_section, "DIPOLE_MOMENTS")
2281 9174 : CALL section_vals_val_get(dipole_section, "DIPOLE_FORM", explicit=explicit)
2282 9174 : IF (explicit) THEN
2283 36 : CALL section_vals_val_get(dipole_section, "DIPOLE_FORM", i_val=t_control%dipole_form)
2284 : ELSE
2285 9138 : t_control%dipole_form = 0
2286 : END IF
2287 9174 : CALL section_vals_val_get(dipole_section, "REFERENCE", i_val=t_control%dipole_reference)
2288 9174 : CALL section_vals_val_get(dipole_section, "REFERENCE_POINT", explicit=explicit)
2289 9174 : IF (explicit) THEN
2290 0 : CALL section_vals_val_get(dipole_section, "REFERENCE_POINT", r_vals=t_control%dipole_ref_point)
2291 : ELSE
2292 9174 : NULLIFY (t_control%dipole_ref_point)
2293 9174 : IF (t_control%dipole_form == tddfpt_dipole_length .AND. t_control%dipole_reference == use_mom_ref_user) THEN
2294 0 : CPABORT("User-defined reference point should be given explicitly")
2295 : END IF
2296 : END IF
2297 :
2298 : !SOC subsection
2299 9174 : soc_section => section_vals_get_subs_vals(t_section, "SOC")
2300 9174 : CALL section_vals_get(soc_section, explicit=explicit)
2301 9174 : IF (explicit) THEN
2302 10 : t_control%do_soc = .TRUE.
2303 : END IF
2304 :
2305 : ! MGRID subsection
2306 9174 : mgrid_section => section_vals_get_subs_vals(t_section, "MGRID")
2307 9174 : CALL section_vals_get(mgrid_section, explicit=t_control%mgrid_is_explicit)
2308 :
2309 9174 : IF (t_control%mgrid_is_explicit) THEN
2310 10 : CALL section_vals_val_get(mgrid_section, "NGRIDS", i_val=t_control%mgrid_ngrids, explicit=explicit)
2311 10 : IF (.NOT. explicit) t_control%mgrid_ngrids = SIZE(qs_control%e_cutoff)
2312 :
2313 10 : CALL section_vals_val_get(mgrid_section, "CUTOFF", r_val=t_control%mgrid_cutoff, explicit=explicit)
2314 10 : IF (.NOT. explicit) t_control%mgrid_cutoff = qs_control%cutoff
2315 :
2316 : CALL section_vals_val_get(mgrid_section, "PROGRESSION_FACTOR", &
2317 10 : r_val=t_control%mgrid_progression_factor, explicit=explicit)
2318 10 : IF (explicit) THEN
2319 0 : IF (t_control%mgrid_progression_factor <= 1.0_dp) THEN
2320 : CALL cp_abort(__LOCATION__, &
2321 0 : "Progression factor should be greater then 1.0 to ensure multi-grid ordering")
2322 : END IF
2323 : ELSE
2324 10 : t_control%mgrid_progression_factor = qs_control%progression_factor
2325 : END IF
2326 :
2327 10 : CALL section_vals_val_get(mgrid_section, "COMMENSURATE", l_val=t_control%mgrid_commensurate_mgrids, explicit=explicit)
2328 10 : IF (.NOT. explicit) t_control%mgrid_commensurate_mgrids = qs_control%commensurate_mgrids
2329 10 : IF (t_control%mgrid_commensurate_mgrids) THEN
2330 0 : IF (explicit) THEN
2331 0 : t_control%mgrid_progression_factor = 4.0_dp
2332 : ELSE
2333 0 : t_control%mgrid_progression_factor = qs_control%progression_factor
2334 : END IF
2335 : END IF
2336 :
2337 10 : CALL section_vals_val_get(mgrid_section, "REL_CUTOFF", r_val=t_control%mgrid_relative_cutoff, explicit=explicit)
2338 10 : IF (.NOT. explicit) t_control%mgrid_relative_cutoff = qs_control%relative_cutoff
2339 :
2340 10 : CALL section_vals_val_get(mgrid_section, "MULTIGRID_SET", l_val=multigrid_set, explicit=explicit)
2341 10 : IF (.NOT. explicit) multigrid_set = .FALSE.
2342 10 : IF (multigrid_set) THEN
2343 0 : CALL section_vals_val_get(mgrid_section, "MULTIGRID_CUTOFF", r_vals=t_control%mgrid_e_cutoff)
2344 : ELSE
2345 10 : NULLIFY (t_control%mgrid_e_cutoff)
2346 : END IF
2347 :
2348 10 : CALL section_vals_val_get(mgrid_section, "REALSPACE", l_val=t_control%mgrid_realspace_mgrids, explicit=explicit)
2349 10 : IF (.NOT. explicit) t_control%mgrid_realspace_mgrids = qs_control%realspace_mgrids
2350 :
2351 : CALL section_vals_val_get(mgrid_section, "SKIP_LOAD_BALANCE_DISTRIBUTED", &
2352 10 : l_val=t_control%mgrid_skip_load_balance, explicit=explicit)
2353 10 : IF (.NOT. explicit) t_control%mgrid_skip_load_balance = qs_control%skip_load_balance_distributed
2354 :
2355 10 : IF (ASSOCIATED(t_control%mgrid_e_cutoff)) THEN
2356 0 : IF (SIZE(t_control%mgrid_e_cutoff) /= t_control%mgrid_ngrids) THEN
2357 0 : CPABORT("Inconsistent values for number of multi-grids")
2358 : END IF
2359 :
2360 : ! sort multi-grids in descending order according to their cutoff values
2361 0 : t_control%mgrid_e_cutoff = -t_control%mgrid_e_cutoff
2362 0 : ALLOCATE (inds(t_control%mgrid_ngrids))
2363 0 : CALL sort(t_control%mgrid_e_cutoff, t_control%mgrid_ngrids, inds)
2364 0 : DEALLOCATE (inds)
2365 0 : t_control%mgrid_e_cutoff = -t_control%mgrid_e_cutoff
2366 : END IF
2367 : END IF
2368 :
2369 : ! expand XC subsection (if given explicitly)
2370 9174 : xc_section => section_vals_get_subs_vals(t_section, "XC")
2371 9174 : xc_func => section_vals_get_subs_vals(xc_section, "XC_FUNCTIONAL")
2372 9174 : CALL section_vals_get(xc_func, explicit=explicit)
2373 9174 : IF (explicit) THEN
2374 298 : CALL xc_functionals_expand(xc_func, xc_section)
2375 : END IF
2376 :
2377 : ! sTDA subsection
2378 9174 : stda_section => section_vals_get_subs_vals(t_section, "STDA")
2379 9174 : IF (t_control%kernel == tddfpt_kernel_stda) THEN
2380 132 : t_control%stda_control%hfx_fraction = 0.0_dp
2381 132 : t_control%stda_control%do_exchange = .TRUE.
2382 132 : t_control%stda_control%eps_td_filter = 1.e-10_dp
2383 132 : t_control%stda_control%mn_alpha = -99.0_dp
2384 132 : t_control%stda_control%mn_beta = -99.0_dp
2385 : ! set default for Ewald method (on/off) dependent on periodicity
2386 236 : SELECT CASE (qs_control%periodicity)
2387 : CASE (0)
2388 104 : t_control%stda_control%do_ewald = .FALSE.
2389 : CASE (1)
2390 0 : t_control%stda_control%do_ewald = .TRUE.
2391 : CASE (2)
2392 0 : t_control%stda_control%do_ewald = .TRUE.
2393 : CASE (3)
2394 28 : t_control%stda_control%do_ewald = .TRUE.
2395 : CASE DEFAULT
2396 132 : CPABORT("Illegal value for periodiciy")
2397 : END SELECT
2398 132 : CALL section_vals_get(stda_section, explicit=explicit)
2399 132 : IF (explicit) THEN
2400 116 : CALL section_vals_val_get(stda_section, "HFX_FRACTION", r_val=hfx, explicit=expl)
2401 116 : IF (expl) t_control%stda_control%hfx_fraction = hfx
2402 116 : CALL section_vals_val_get(stda_section, "EPS_TD_FILTER", r_val=filter, explicit=expl)
2403 116 : IF (expl) t_control%stda_control%eps_td_filter = filter
2404 116 : CALL section_vals_val_get(stda_section, "DO_EWALD", l_val=do_ewald, explicit=expl)
2405 116 : IF (expl) t_control%stda_control%do_ewald = do_ewald
2406 116 : CALL section_vals_val_get(stda_section, "DO_EXCHANGE", l_val=do_exchange, explicit=expl)
2407 116 : IF (expl) t_control%stda_control%do_exchange = do_exchange
2408 116 : CALL section_vals_val_get(stda_section, "MATAGA_NISHIMOTO_CEXP", r_val=fval)
2409 116 : t_control%stda_control%mn_alpha = fval
2410 116 : CALL section_vals_val_get(stda_section, "MATAGA_NISHIMOTO_XEXP", r_val=fval)
2411 116 : t_control%stda_control%mn_beta = fval
2412 : END IF
2413 132 : CALL section_vals_val_get(stda_section, "COULOMB_SR_CUT", r_val=fval)
2414 132 : t_control%stda_control%coulomb_sr_cut = fval
2415 132 : CALL section_vals_val_get(stda_section, "COULOMB_SR_EPS", r_val=fval)
2416 132 : t_control%stda_control%coulomb_sr_eps = fval
2417 : END IF
2418 :
2419 9174 : CALL timestop(handle)
2420 9174 : END SUBROUTINE read_tddfpt2_control
2421 :
2422 : ! **************************************************************************************************
2423 : !> \brief Write the DFT control parameters to the output unit.
2424 : !> \param dft_control ...
2425 : !> \param dft_section ...
2426 : ! **************************************************************************************************
2427 15718 : SUBROUTINE write_dft_control(dft_control, dft_section)
2428 : TYPE(dft_control_type), POINTER :: dft_control
2429 : TYPE(section_vals_type), POINTER :: dft_section
2430 :
2431 : CHARACTER(len=*), PARAMETER :: routineN = 'write_dft_control'
2432 :
2433 : CHARACTER(LEN=20) :: tmpStr
2434 : INTEGER :: handle, i, i_rep, max_mtlr_iter, n_rep, &
2435 : output_unit
2436 : REAL(kind=dp) :: density_cut, density_smooth_cut_range, &
2437 : eps_u_j_loop, gradient_cut, tau_cut
2438 : TYPE(cp_logger_type), POINTER :: logger
2439 : TYPE(enumeration_type), POINTER :: enum
2440 : TYPE(keyword_type), POINTER :: keyword
2441 : TYPE(section_type), POINTER :: section
2442 : TYPE(section_vals_type), POINTER :: xc_section
2443 :
2444 10654 : IF (dft_control%qs_control%semi_empirical) RETURN
2445 8124 : IF (dft_control%qs_control%dftb) RETURN
2446 7826 : IF (dft_control%qs_control%xtb) THEN
2447 1232 : CALL write_xtb_control(dft_control%qs_control%xtb_control, dft_section)
2448 1232 : RETURN
2449 : END IF
2450 6594 : CALL timeset(routineN, handle)
2451 :
2452 6594 : NULLIFY (logger)
2453 6594 : logger => cp_get_default_logger()
2454 :
2455 : output_unit = cp_print_key_unit_nr(logger, dft_section, &
2456 6594 : "PRINT%DFT_CONTROL_PARAMETERS", extension=".Log")
2457 :
2458 6594 : IF (output_unit > 0) THEN
2459 :
2460 1557 : xc_section => section_vals_get_subs_vals(dft_section, "XC")
2461 :
2462 1557 : IF (dft_control%uks) THEN
2463 : WRITE (UNIT=output_unit, FMT="(/,T2,A,T78,A)") &
2464 442 : "DFT| Spin unrestricted (spin-polarized) Kohn-Sham calculation", "UKS"
2465 1115 : ELSE IF (dft_control%roks) THEN
2466 : WRITE (UNIT=output_unit, FMT="(/,T2,A,T77,A)") &
2467 15 : "DFT| Spin restricted open Kohn-Sham calculation", "ROKS"
2468 : ELSE
2469 : WRITE (UNIT=output_unit, FMT="(/,T2,A,T78,A)") &
2470 1100 : "DFT| Spin restricted Kohn-Sham (RKS) calculation", "RKS"
2471 : END IF
2472 :
2473 : WRITE (UNIT=output_unit, FMT="(T2,A,T76,I5)") &
2474 1557 : "DFT| Multiplicity", dft_control%multiplicity
2475 : WRITE (UNIT=output_unit, FMT="(T2,A,T76,I5)") &
2476 1557 : "DFT| Number of spin states", dft_control%nspins
2477 :
2478 : WRITE (UNIT=output_unit, FMT="(T2,A,T76,I5)") &
2479 1557 : "DFT| Charge", dft_control%charge
2480 :
2481 1557 : IF (dft_control%sic_method_id /= sic_none) CALL cite_reference(VandeVondele2005b)
2482 3100 : SELECT CASE (dft_control%sic_method_id)
2483 : CASE (sic_none)
2484 1543 : tmpstr = "NO"
2485 : CASE (sic_mauri_spz)
2486 6 : tmpstr = "SPZ/MAURI SIC"
2487 : CASE (sic_mauri_us)
2488 3 : tmpstr = "US/MAURI SIC"
2489 : CASE (sic_ad)
2490 3 : tmpstr = "AD SIC"
2491 : CASE (sic_eo)
2492 2 : tmpstr = "Explicit Orbital SIC"
2493 : CASE DEFAULT
2494 : ! fix throughout the cp2k for this option
2495 1557 : CPABORT("SIC option unknown")
2496 : END SELECT
2497 :
2498 : WRITE (UNIT=output_unit, FMT="(T2,A,T61,A20)") &
2499 1557 : "DFT| Self-interaction correction (SIC)", ADJUSTR(TRIM(tmpstr))
2500 :
2501 1557 : IF (dft_control%sic_method_id /= sic_none) THEN
2502 : WRITE (UNIT=output_unit, FMT="(T2,A,T66,ES15.6)") &
2503 14 : "DFT| SIC scaling parameter a", dft_control%sic_scaling_a, &
2504 28 : "DFT| SIC scaling parameter b", dft_control%sic_scaling_b
2505 : END IF
2506 :
2507 1557 : IF (dft_control%sic_method_id == sic_eo) THEN
2508 2 : IF (dft_control%sic_list_id == sic_list_all) THEN
2509 : WRITE (UNIT=output_unit, FMT="(T2,A,T66,A)") &
2510 1 : "DFT| SIC orbitals", "ALL"
2511 : END IF
2512 2 : IF (dft_control%sic_list_id == sic_list_unpaired) THEN
2513 : WRITE (UNIT=output_unit, FMT="(T2,A,T66,A)") &
2514 1 : "DFT| SIC orbitals", "UNPAIRED"
2515 : END IF
2516 : END IF
2517 :
2518 1557 : CALL section_vals_val_get(xc_section, "density_cutoff", r_val=density_cut)
2519 1557 : CALL section_vals_val_get(xc_section, "gradient_cutoff", r_val=gradient_cut)
2520 1557 : CALL section_vals_val_get(xc_section, "tau_cutoff", r_val=tau_cut)
2521 1557 : CALL section_vals_val_get(xc_section, "density_smooth_cutoff_range", r_val=density_smooth_cut_range)
2522 :
2523 : WRITE (UNIT=output_unit, FMT="(T2,A,T66,ES15.6)") &
2524 1557 : "DFT| Cutoffs: density ", density_cut, &
2525 1557 : "DFT| gradient", gradient_cut, &
2526 1557 : "DFT| tau ", tau_cut, &
2527 3114 : "DFT| cutoff_smoothing_range", density_smooth_cut_range
2528 : CALL section_vals_val_get(xc_section, "XC_GRID%XC_SMOOTH_RHO", &
2529 1557 : c_val=tmpStr)
2530 : WRITE (output_unit, '( A, T61, A )') &
2531 1557 : " DFT| XC density smoothing ", ADJUSTR(tmpStr)
2532 : CALL section_vals_val_get(xc_section, "XC_GRID%XC_DERIV", &
2533 1557 : c_val=tmpStr)
2534 : WRITE (output_unit, '( A, T61, A )') &
2535 1557 : " DFT| XC derivatives ", ADJUSTR(tmpStr)
2536 1557 : IF (dft_control%dft_plus_u) THEN
2537 19 : NULLIFY (enum, keyword, section)
2538 19 : CALL create_dft_section(section)
2539 19 : keyword => section_get_keyword(section, "PLUS_U_METHOD")
2540 19 : CALL keyword_get(keyword, enum=enum)
2541 : WRITE (UNIT=output_unit, FMT="(/,T2,A,T41,A40)") &
2542 19 : "DFT+U| Method", ADJUSTR(TRIM(enum_i2c(enum, dft_control%plus_u_method_id)))
2543 : WRITE (UNIT=output_unit, FMT="(T2,A)") &
2544 19 : "DFT+U| Check atomic kind information for details"
2545 19 : IF (dft_control%mtlr_u_j) THEN
2546 0 : CALL section_vals_val_get(dft_section, "EPS_U_J_LOOP", r_val=eps_u_j_loop)
2547 : WRITE (UNIT=output_unit, FMT="(T2,A,T67,ES14.7E3)") &
2548 0 : "MTLR U J| EPS_U_J_LOOP", eps_u_j_loop
2549 0 : CALL section_vals_val_get(dft_section, "MAX_MTLR_LOOP", i_val=max_mtlr_iter)
2550 : WRITE (UNIT=output_unit, FMT="(T2,A,T67,I10)") &
2551 0 : "MTLR U J| MAX_MTLR_LOOP", max_mtlr_iter
2552 : END IF
2553 19 : CALL section_release(section)
2554 : END IF
2555 :
2556 1557 : WRITE (UNIT=output_unit, FMT="(A)") ""
2557 1557 : CALL xc_write(output_unit, xc_section, dft_control%lsd)
2558 :
2559 1557 : IF (dft_control%apply_period_efield) THEN
2560 6 : WRITE (UNIT=output_unit, FMT="(A)") ""
2561 6 : IF (dft_control%period_efield%displacement_field) THEN
2562 : WRITE (UNIT=output_unit, FMT="(T2,A)") &
2563 0 : "PERIODIC_EFIELD| Use displacement field formulation"
2564 : WRITE (UNIT=output_unit, FMT="(T2,A,T66,1X,ES14.6)") &
2565 0 : "PERIODIC_EFIELD| Displacement field filter: x", &
2566 0 : dft_control%period_efield%d_filter(1), &
2567 0 : "PERIODIC_EFIELD| y", &
2568 0 : dft_control%period_efield%d_filter(2), &
2569 0 : "PERIODIC_EFIELD| z", &
2570 0 : dft_control%period_efield%d_filter(3)
2571 : END IF
2572 : WRITE (UNIT=output_unit, FMT="(T2,A,T66,1X,ES14.6)") &
2573 6 : "PERIODIC_EFIELD| Polarisation vector: x", &
2574 6 : dft_control%period_efield%polarisation(1), &
2575 6 : "PERIODIC_EFIELD| y", &
2576 6 : dft_control%period_efield%polarisation(2), &
2577 6 : "PERIODIC_EFIELD| z", &
2578 12 : dft_control%period_efield%polarisation(3)
2579 :
2580 : WRITE (UNIT=output_unit, FMT="(T2,A,T66,1X,I14)") &
2581 6 : "PERIODIC_EFIELD| Start Frame:", &
2582 6 : dft_control%period_efield%start_frame, &
2583 6 : "PERIODIC_EFIELD| End Frame:", &
2584 12 : dft_control%period_efield%end_frame
2585 :
2586 6 : IF (ALLOCATED(dft_control%period_efield%strength_list)) THEN
2587 : WRITE (UNIT=output_unit, FMT="(T2,A,T66,1X,I14)") &
2588 2 : "PERIODIC_EFIELD| Number of Intensities:", &
2589 4 : SIZE(dft_control%period_efield%strength_list)
2590 : WRITE (UNIT=output_unit, FMT="(T2,A,I10,T66,1X,ES14.6)") &
2591 2 : "PERIODIC_EFIELD| Intensity List [a.u.] ", &
2592 4 : 1, dft_control%period_efield%strength_list(1)
2593 24 : DO i = 2, SIZE(dft_control%period_efield%strength_list)
2594 : WRITE (UNIT=output_unit, FMT="(T2,A,I10,T66,1X,ES14.6)") &
2595 22 : "PERIODIC_EFIELD| ", &
2596 46 : i, dft_control%period_efield%strength_list(i)
2597 : END DO
2598 : ELSE
2599 : WRITE (UNIT=output_unit, FMT="(T2,A,T66,1X,ES14.6)") &
2600 4 : "PERIODIC_EFIELD| Intensity [a.u.]:", &
2601 8 : dft_control%period_efield%strength
2602 : END IF
2603 :
2604 24 : IF (NORM2(dft_control%period_efield%polarisation) < EPSILON(0.0_dp)) THEN
2605 0 : CPABORT("Invalid (too small) polarisation vector specified for PERIODIC_EFIELD")
2606 : END IF
2607 : END IF
2608 :
2609 1557 : IF (dft_control%do_sccs) THEN
2610 : WRITE (UNIT=output_unit, FMT="(/,T2,A)") &
2611 5 : "SCCS| Self-consistent continuum solvation model"
2612 : WRITE (UNIT=output_unit, FMT="(T2,A,T61,ES20.6)") &
2613 5 : "SCCS| Relative permittivity of the solvent (medium)", &
2614 5 : dft_control%sccs_control%epsilon_solvent, &
2615 5 : "SCCS| Absolute permittivity [a.u.]", &
2616 10 : dft_control%sccs_control%epsilon_solvent/fourpi
2617 9 : SELECT CASE (dft_control%sccs_control%method_id)
2618 : CASE (sccs_andreussi)
2619 : WRITE (UNIT=output_unit, FMT="(T2,A,/,(T2,A,T61,ES20.6))") &
2620 4 : "SCCS| Dielectric function proposed by Andreussi et al.", &
2621 4 : "SCCS| rho_max", dft_control%sccs_control%rho_max, &
2622 8 : "SCCS| rho_min", dft_control%sccs_control%rho_min
2623 : CASE (sccs_fattebert_gygi)
2624 : WRITE (UNIT=output_unit, FMT="(T2,A,/,(T2,A,T61,ES20.6))") &
2625 1 : "SCCS| Dielectric function proposed by Fattebert and Gygi", &
2626 1 : "SCCS| beta", dft_control%sccs_control%beta, &
2627 2 : "SCCS| rho_zero", dft_control%sccs_control%rho_zero
2628 : CASE (sccs_saa_andreussi)
2629 : WRITE (UNIT=output_unit, FMT="(T2,A,/,A,/,(T2,A,T61,ES20.6))") &
2630 0 : "SCCS| Dielectric function of the solvent aware algorithm", &
2631 0 : "SCCS| proposed by Andreussi et al.", &
2632 0 : "SCCS| rho_max", dft_control%sccs_control%rho_max, &
2633 0 : "SCCS| rho_min", dft_control%sccs_control%rho_min, &
2634 0 : "SCCS| f0", dft_control%sccs_control%f0, &
2635 0 : "SCCS| delta_eta", dft_control%sccs_control%delta_eta, &
2636 0 : "SCCS| alpha_zeta", dft_control%sccs_control%alpha_zeta, &
2637 0 : "SCCS| delta_zeta", dft_control%sccs_control%delta_zeta, &
2638 0 : "SCCS| R_solv", dft_control%sccs_control%R_solv
2639 : CASE DEFAULT
2640 5 : CPABORT("Invalid SCCS model specified. Please, check your input!")
2641 : END SELECT
2642 6 : SELECT CASE (dft_control%sccs_control%derivative_method)
2643 : CASE (sccs_derivative_fft)
2644 : WRITE (UNIT=output_unit, FMT="(T2,A,T46,A35)") &
2645 1 : "SCCS| Numerical derivative calculation", &
2646 2 : ADJUSTR("FFT")
2647 : CASE (sccs_derivative_cd3)
2648 : WRITE (UNIT=output_unit, FMT="(T2,A,T46,A35)") &
2649 0 : "SCCS| Numerical derivative calculation", &
2650 0 : ADJUSTR("3-point stencil central differences")
2651 : CASE (sccs_derivative_cd5)
2652 : WRITE (UNIT=output_unit, FMT="(T2,A,T46,A35)") &
2653 4 : "SCCS| Numerical derivative calculation", &
2654 8 : ADJUSTR("5-point stencil central differences")
2655 : CASE (sccs_derivative_cd7)
2656 : WRITE (UNIT=output_unit, FMT="(T2,A,T46,A35)") &
2657 0 : "SCCS| Numerical derivative calculation", &
2658 0 : ADJUSTR("7-point stencil central differences")
2659 : CASE DEFAULT
2660 : CALL cp_abort(__LOCATION__, &
2661 : "Invalid derivative method specified for SCCS model. "// &
2662 5 : "Please, check your input!")
2663 : END SELECT
2664 : WRITE (UNIT=output_unit, FMT="(T2,A,T61,ES20.6)") &
2665 5 : "SCCS| Repulsion parameter alpha [mN/m] = [dyn/cm]", &
2666 10 : cp_unit_from_cp2k(dft_control%sccs_control%alpha_solvent, "mN/m")
2667 : WRITE (UNIT=output_unit, FMT="(T2,A,T61,ES20.6)") &
2668 5 : "SCCS| Dispersion parameter beta [GPa]", &
2669 10 : cp_unit_from_cp2k(dft_control%sccs_control%beta_solvent, "GPa")
2670 : WRITE (UNIT=output_unit, FMT="(T2,A,T61,ES20.6)") &
2671 5 : "SCCS| Surface tension gamma [mN/m] = [dyn/cm]", &
2672 10 : cp_unit_from_cp2k(dft_control%sccs_control%gamma_solvent, "mN/m")
2673 : WRITE (UNIT=output_unit, FMT="(T2,A,T61,ES20.6)") &
2674 5 : "SCCS| Mixing parameter applied during the iteration cycle", &
2675 10 : dft_control%sccs_control%mixing
2676 : WRITE (UNIT=output_unit, FMT="(T2,A,T61,ES20.6)") &
2677 5 : "SCCS| Tolerance for the convergence of the SCCS iteration cycle", &
2678 10 : dft_control%sccs_control%eps_sccs
2679 : WRITE (UNIT=output_unit, FMT="(T2,A,T61,I20)") &
2680 5 : "SCCS| Maximum number of iteration steps", &
2681 10 : dft_control%sccs_control%max_iter
2682 : WRITE (UNIT=output_unit, FMT="(T2,A,T61,ES20.6)") &
2683 5 : "SCCS| SCF convergence threshold for starting the SCCS iteration", &
2684 10 : dft_control%sccs_control%eps_scf
2685 : WRITE (UNIT=output_unit, FMT="(T2,A,T61,ES20.6)") &
2686 5 : "SCCS| Numerical increment for the cavity surface calculation", &
2687 10 : dft_control%sccs_control%delta_rho
2688 : END IF
2689 :
2690 1557 : WRITE (UNIT=output_unit, FMT="(A)") ""
2691 :
2692 : END IF
2693 :
2694 6594 : IF (dft_control%hairy_probes .EQV. .TRUE.) THEN
2695 4 : n_rep = SIZE(dft_control%probe)
2696 4 : IF (output_unit > 0) THEN
2697 6 : DO i_rep = 1, n_rep
2698 : WRITE (UNIT=output_unit, FMT="(T2,A,I5)") &
2699 4 : "HP | hair probe set", i_rep
2700 : WRITE (UNIT=output_unit, FMT="(T2,A,T61,*(I5))") &
2701 4 : "HP| atom indexes", &
2702 12 : (dft_control%probe(i_rep)%atom_ids(i), i=1, dft_control%probe(i_rep)%natoms)
2703 : WRITE (UNIT=output_unit, FMT="(T2,A,T61,ES20.6)") &
2704 4 : "HP| potential", dft_control%probe(i_rep)%mu
2705 : WRITE (UNIT=output_unit, FMT="(T2,A,T61,F20.2)") &
2706 4 : "HP| temperature", dft_control%probe(i_rep)%T
2707 : WRITE (UNIT=output_unit, FMT="(T2,A,T61,ES20.6)") &
2708 6 : "HP| eps_hp", dft_control%probe(i_rep)%eps_hp
2709 : END DO
2710 : END IF
2711 : END IF
2712 :
2713 : CALL cp_print_key_finished_output(output_unit, logger, dft_section, &
2714 6594 : "PRINT%DFT_CONTROL_PARAMETERS")
2715 :
2716 6594 : CALL timestop(handle)
2717 :
2718 : END SUBROUTINE write_dft_control
2719 :
2720 : ! **************************************************************************************************
2721 : !> \brief Write the ADMM control parameters to the output unit.
2722 : !> \param admm_control ...
2723 : !> \param dft_section ...
2724 : ! **************************************************************************************************
2725 524 : SUBROUTINE write_admm_control(admm_control, dft_section)
2726 : TYPE(admm_control_type), POINTER :: admm_control
2727 : TYPE(section_vals_type), POINTER :: dft_section
2728 :
2729 : INTEGER :: iounit
2730 : TYPE(cp_logger_type), POINTER :: logger
2731 :
2732 524 : NULLIFY (logger)
2733 524 : logger => cp_get_default_logger()
2734 :
2735 : iounit = cp_print_key_unit_nr(logger, dft_section, &
2736 524 : "PRINT%DFT_CONTROL_PARAMETERS", extension=".Log")
2737 :
2738 524 : IF (iounit > 0) THEN
2739 :
2740 257 : SELECT CASE (admm_control%admm_type)
2741 : CASE (no_admm_type)
2742 124 : WRITE (UNIT=iounit, FMT="(/,T2,A,T77,A)") "ADMM| Specific ADMM type specified", "NONE"
2743 : CASE (admm1_type)
2744 2 : WRITE (UNIT=iounit, FMT="(/,T2,A,T76,A)") "ADMM| Specific ADMM type specified", "ADMM1"
2745 : CASE (admm2_type)
2746 1 : WRITE (UNIT=iounit, FMT="(/,T2,A,T76,A)") "ADMM| Specific ADMM type specified", "ADMM2"
2747 : CASE (admms_type)
2748 4 : WRITE (UNIT=iounit, FMT="(/,T2,A,T76,A)") "ADMM| Specific ADMM type specified", "ADMMS"
2749 : CASE (admmp_type)
2750 1 : WRITE (UNIT=iounit, FMT="(/,T2,A,T76,A)") "ADMM| Specific ADMM type specified", "ADMMP"
2751 : CASE (admmq_type)
2752 1 : WRITE (UNIT=iounit, FMT="(/,T2,A,T76,A)") "ADMM| Specific ADMM type specified", "ADMMQ"
2753 : CASE DEFAULT
2754 133 : CPABORT("admm_type")
2755 : END SELECT
2756 :
2757 217 : SELECT CASE (admm_control%purification_method)
2758 : CASE (do_admm_purify_none)
2759 84 : WRITE (UNIT=iounit, FMT="(T2,A,T77,A)") "ADMM| Density matrix purification method", "NONE"
2760 : CASE (do_admm_purify_cauchy)
2761 9 : WRITE (UNIT=iounit, FMT="(T2,A,T75,A)") "ADMM| Density matrix purification method", "Cauchy"
2762 : CASE (do_admm_purify_cauchy_subspace)
2763 5 : WRITE (UNIT=iounit, FMT="(T2,A,T66,A)") "ADMM| Density matrix purification method", "Cauchy subspace"
2764 : CASE (do_admm_purify_mo_diag)
2765 25 : WRITE (UNIT=iounit, FMT="(T2,A,T63,A)") "ADMM| Density matrix purification method", "MO diagonalization"
2766 : CASE (do_admm_purify_mo_no_diag)
2767 3 : WRITE (UNIT=iounit, FMT="(T2,A,T71,A)") "ADMM| Density matrix purification method", "MO no diag"
2768 : CASE (do_admm_purify_mcweeny)
2769 1 : WRITE (UNIT=iounit, FMT="(T2,A,T74,A)") "ADMM| Density matrix purification method", "McWeeny"
2770 : CASE (do_admm_purify_none_dm)
2771 6 : WRITE (UNIT=iounit, FMT="(T2,A,T73,A)") "ADMM| Density matrix purification method", "NONE(DM)"
2772 : CASE DEFAULT
2773 133 : CPABORT("admm_purification_method")
2774 : END SELECT
2775 :
2776 238 : SELECT CASE (admm_control%method)
2777 : CASE (do_admm_basis_projection)
2778 105 : WRITE (UNIT=iounit, FMT="(T2,A)") "ADMM| Orbital projection on ADMM basis"
2779 : CASE (do_admm_blocking_purify_full)
2780 3 : WRITE (UNIT=iounit, FMT="(T2,A)") "ADMM| Blocked Fock matrix projection with full purification"
2781 : CASE (do_admm_blocked_projection)
2782 6 : WRITE (UNIT=iounit, FMT="(T2,A)") "ADMM| Blocked Fock matrix projection"
2783 : CASE (do_admm_charge_constrained_projection)
2784 19 : WRITE (UNIT=iounit, FMT="(T2,A)") "ADMM| Orbital projection with charge constrain"
2785 : CASE DEFAULT
2786 133 : CPABORT("admm method")
2787 : END SELECT
2788 :
2789 154 : SELECT CASE (admm_control%scaling_model)
2790 : CASE (do_admm_exch_scaling_none)
2791 : CASE (do_admm_exch_scaling_merlot)
2792 21 : WRITE (UNIT=iounit, FMT="(T2,A)") "ADMM| Use Merlot (2014) scaling model"
2793 : CASE DEFAULT
2794 133 : CPABORT("admm scaling_model")
2795 : END SELECT
2796 :
2797 133 : WRITE (UNIT=iounit, FMT="(T2,A,T61,G20.10)") "ADMM| eps_filter", admm_control%eps_filter
2798 :
2799 145 : SELECT CASE (admm_control%aux_exch_func)
2800 : CASE (do_admm_aux_exch_func_none)
2801 12 : WRITE (UNIT=iounit, FMT="(T2,A)") "ADMM| No exchange functional correction term used"
2802 : CASE (do_admm_aux_exch_func_default, do_admm_aux_exch_func_default_libxc)
2803 90 : WRITE (UNIT=iounit, FMT="(T2,A,T74,A)") "ADMM| Exchange functional in correction term", "(W)PBEX"
2804 : CASE (do_admm_aux_exch_func_pbex, do_admm_aux_exch_func_pbex_libxc)
2805 22 : WRITE (UNIT=iounit, FMT="(T2,A,T77,A)") "ADMM| Exchange functional in correction term", "PBEX"
2806 : CASE (do_admm_aux_exch_func_opt, do_admm_aux_exch_func_opt_libxc)
2807 8 : WRITE (UNIT=iounit, FMT="(T2,A,T77,A)") "ADMM| Exchange functional in correction term", "OPTX"
2808 : CASE (do_admm_aux_exch_func_bee, do_admm_aux_exch_func_bee_libxc)
2809 1 : WRITE (UNIT=iounit, FMT="(T2,A,T74,A)") "ADMM| Exchange functional in correction term", "Becke88"
2810 : CASE (do_admm_aux_exch_func_sx_libxc)
2811 0 : WRITE (UNIT=iounit, FMT="(T2,A,T74,A)") "ADMM| Exchange functional in correction term", "SlaterX"
2812 : CASE DEFAULT
2813 133 : CPABORT("admm aux_exch_func")
2814 : END SELECT
2815 :
2816 133 : WRITE (UNIT=iounit, FMT="(A)") ""
2817 :
2818 : END IF
2819 :
2820 : CALL cp_print_key_finished_output(iounit, logger, dft_section, &
2821 524 : "PRINT%DFT_CONTROL_PARAMETERS")
2822 524 : END SUBROUTINE write_admm_control
2823 :
2824 : ! **************************************************************************************************
2825 : !> \brief Write the xTB control parameters to the output unit.
2826 : !> \param xtb_control ...
2827 : !> \param dft_section ...
2828 : ! **************************************************************************************************
2829 1232 : SUBROUTINE write_xtb_control(xtb_control, dft_section)
2830 : TYPE(xtb_control_type), POINTER :: xtb_control
2831 : TYPE(section_vals_type), POINTER :: dft_section
2832 :
2833 : CHARACTER(len=*), PARAMETER :: routineN = 'write_xtb_control'
2834 :
2835 : CHARACTER(LEN=16) :: scc_mixer_name, solver_name
2836 : INTEGER :: handle, output_unit
2837 : TYPE(cp_logger_type), POINTER :: logger
2838 :
2839 1232 : CALL timeset(routineN, handle)
2840 1232 : NULLIFY (logger)
2841 1232 : logger => cp_get_default_logger()
2842 :
2843 : output_unit = cp_print_key_unit_nr(logger, dft_section, &
2844 1232 : "PRINT%DFT_CONTROL_PARAMETERS", extension=".Log")
2845 :
2846 1232 : IF (output_unit > 0) THEN
2847 :
2848 : WRITE (UNIT=output_unit, FMT="(/,T2,A,T31,A50)") &
2849 121 : "xTB| Parameter file", ADJUSTR(TRIM(xtb_control%parameter_file_name))
2850 : WRITE (UNIT=output_unit, FMT="(T2,A,T71,I10)") &
2851 121 : "xTB| Basis expansion STO-NG", xtb_control%sto_ng
2852 : WRITE (UNIT=output_unit, FMT="(T2,A,T71,I10)") &
2853 121 : "xTB| Basis expansion STO-NG for Hydrogen", xtb_control%h_sto_ng
2854 : WRITE (UNIT=output_unit, FMT="(T2,A,T71,E10.4)") &
2855 121 : "xTB| Repulsive pair potential accuracy", xtb_control%eps_pair
2856 : WRITE (UNIT=output_unit, FMT="(T2,A,T71,F10.6)") &
2857 121 : "xTB| Repulsive enhancement factor", xtb_control%enscale
2858 : WRITE (UNIT=output_unit, FMT="(T2,A,T71,L10)") &
2859 121 : "xTB| Halogen interaction potential", xtb_control%xb_interaction
2860 : WRITE (UNIT=output_unit, FMT="(T2,A,T71,F10.3)") &
2861 121 : "xTB| Halogen interaction potential cutoff radius", xtb_control%xb_radius
2862 : WRITE (UNIT=output_unit, FMT="(T2,A,T71,L10)") &
2863 121 : "xTB| Nonbonded interactions", xtb_control%do_nonbonded
2864 121 : SELECT CASE (xtb_control%vdw_type)
2865 : CASE (xtb_vdw_type_none)
2866 0 : WRITE (UNIT=output_unit, FMT="(T2,A)") "xTB| No vdW potential selected"
2867 : CASE (xtb_vdw_type_d3)
2868 120 : WRITE (UNIT=output_unit, FMT="(T2,A,T72,A)") "xTB| vdW potential type:", "DFTD3(BJ)"
2869 : WRITE (UNIT=output_unit, FMT="(T2,A,T31,A50)") &
2870 120 : "xTB| D3 Dispersion: Parameter file", ADJUSTR(TRIM(xtb_control%dispersion_parameter_file))
2871 : CASE (xtb_vdw_type_d4)
2872 1 : WRITE (UNIT=output_unit, FMT="(T2,A,T76,A)") "xTB| vdW potential type:", "DFTD4"
2873 : WRITE (UNIT=output_unit, FMT="(T2,A,T31,A50)") &
2874 1 : "xTB| D4 Dispersion: Parameter file", ADJUSTR(TRIM(xtb_control%dispersion_parameter_file))
2875 : CASE DEFAULT
2876 121 : CPABORT("vdw type")
2877 : END SELECT
2878 : WRITE (UNIT=output_unit, FMT="(T2,A,T51,3F10.3)") &
2879 121 : "xTB| Huckel constants ks kp kd", xtb_control%ks, xtb_control%kp, xtb_control%kd
2880 : WRITE (UNIT=output_unit, FMT="(T2,A,T61,2F10.3)") &
2881 121 : "xTB| Huckel constants ksp k2sh", xtb_control%ksp, xtb_control%k2sh
2882 : WRITE (UNIT=output_unit, FMT="(T2,A,T71,F10.3)") &
2883 121 : "xTB| Mataga-Nishimoto exponent", xtb_control%kg
2884 : WRITE (UNIT=output_unit, FMT="(T2,A,T71,F10.3)") &
2885 121 : "xTB| Repulsion potential exponent", xtb_control%kf
2886 : WRITE (UNIT=output_unit, FMT="(T2,A,T51,3F10.3)") &
2887 121 : "xTB| Coordination number scaling kcn(s) kcn(p) kcn(d)", &
2888 242 : xtb_control%kcns, xtb_control%kcnp, xtb_control%kcnd
2889 : WRITE (UNIT=output_unit, FMT="(T2,A,T71,F10.3)") &
2890 121 : "xTB| Electronegativity scaling", xtb_control%ken
2891 : WRITE (UNIT=output_unit, FMT="(T2,A,T61,2F10.3)") &
2892 121 : "xTB| Halogen potential scaling kxr kx2", xtb_control%kxr, xtb_control%kx2
2893 210 : SELECT CASE (xtb_control%tblite_scc_mixer)
2894 : CASE (tblite_scc_mixer_auto)
2895 89 : scc_mixer_name = "AUTO"
2896 : CASE (tblite_scc_mixer_tblite)
2897 10 : scc_mixer_name = "TBLITE"
2898 : CASE (tblite_scc_mixer_cp2k)
2899 4 : scc_mixer_name = "CP2K"
2900 : CASE (tblite_scc_mixer_none)
2901 18 : scc_mixer_name = "NONE"
2902 : CASE DEFAULT
2903 121 : CPABORT("Unknown tblite SCC mixer")
2904 : END SELECT
2905 242 : SELECT CASE (xtb_control%tblite_mixer_solver)
2906 : CASE (tblite_solver_gvd)
2907 121 : solver_name = "GVD"
2908 : CASE (tblite_solver_gvr)
2909 0 : solver_name = "GVR"
2910 : CASE DEFAULT
2911 121 : CPABORT("Unknown tblite SCC mixer solver")
2912 : END SELECT
2913 : WRITE (UNIT=output_unit, FMT="(T2,A,T72,A)") &
2914 121 : "xTB| SCC mixer:", TRIM(scc_mixer_name)
2915 : WRITE (UNIT=output_unit, FMT="(T2,A,T71,ES10.3)") &
2916 121 : "xTB| tblite accuracy:", xtb_control%tblite_accuracy
2917 121 : IF (LEN_TRIM(xtb_control%tblite_param_file) > 0) THEN
2918 : WRITE (UNIT=output_unit, FMT="(T2,A,T33,A)") &
2919 0 : "xTB| tblite parameter file:", TRIM(xtb_control%tblite_param_file)
2920 : END IF
2921 : WRITE (UNIT=output_unit, FMT="(T2,A,T71,F10.3)") &
2922 121 : "xTB| tblite SCC mixer damping:", xtb_control%tblite_mixer_damping
2923 : WRITE (UNIT=output_unit, FMT="(T2,A,T71,I10)") &
2924 121 : "xTB| tblite SCC mixer iterations:", xtb_control%tblite_mixer_iterations
2925 : WRITE (UNIT=output_unit, FMT="(T2,A,T71,I10)") &
2926 121 : "xTB| tblite SCC mixer memory:", xtb_control%tblite_mixer_memory
2927 : WRITE (UNIT=output_unit, FMT="(T2,A,T72,A)") &
2928 121 : "xTB| tblite SCC mixer solver:", TRIM(solver_name)
2929 : WRITE (UNIT=output_unit, FMT="(T2,A,T71,ES10.3)") &
2930 121 : "xTB| tblite SCC mixer omega0:", xtb_control%tblite_mixer_omega0
2931 : WRITE (UNIT=output_unit, FMT="(T2,A,T71,ES10.3)") &
2932 121 : "xTB| tblite SCC mixer min weight:", xtb_control%tblite_mixer_min_weight
2933 : WRITE (UNIT=output_unit, FMT="(T2,A,T71,ES10.3)") &
2934 121 : "xTB| tblite SCC mixer max weight:", xtb_control%tblite_mixer_max_weight
2935 : WRITE (UNIT=output_unit, FMT="(T2,A,T71,ES10.3)") &
2936 121 : "xTB| tblite SCC mixer weight factor:", xtb_control%tblite_mixer_weight_factor
2937 121 : WRITE (UNIT=output_unit, FMT="(/)")
2938 :
2939 : END IF
2940 :
2941 : CALL cp_print_key_finished_output(output_unit, logger, dft_section, &
2942 1232 : "PRINT%DFT_CONTROL_PARAMETERS")
2943 :
2944 1232 : CALL timestop(handle)
2945 :
2946 1232 : END SUBROUTINE write_xtb_control
2947 :
2948 : ! **************************************************************************************************
2949 : !> \brief Purpose: Write the QS control parameters to the output unit.
2950 : !> \param qs_control ...
2951 : !> \param dft_section ...
2952 : ! **************************************************************************************************
2953 15718 : SUBROUTINE write_qs_control(qs_control, dft_section)
2954 : TYPE(qs_control_type), INTENT(IN) :: qs_control
2955 : TYPE(section_vals_type), POINTER :: dft_section
2956 :
2957 : CHARACTER(len=*), PARAMETER :: routineN = 'write_qs_control'
2958 :
2959 : CHARACTER(len=20) :: method, quadrature
2960 : INTEGER :: handle, i, igrid_level, ngrid_level, &
2961 : output_unit
2962 : TYPE(cp_logger_type), POINTER :: logger
2963 : TYPE(ddapc_restraint_type), POINTER :: ddapc_restraint_control
2964 : TYPE(enumeration_type), POINTER :: enum
2965 : TYPE(keyword_type), POINTER :: keyword
2966 : TYPE(section_type), POINTER :: qs_section
2967 : TYPE(section_vals_type), POINTER :: print_section_vals, qs_section_vals
2968 :
2969 10654 : IF (qs_control%semi_empirical) RETURN
2970 8124 : IF (qs_control%dftb) RETURN
2971 7826 : IF (qs_control%xtb) RETURN
2972 6594 : CALL timeset(routineN, handle)
2973 6594 : NULLIFY (logger, print_section_vals, qs_section, qs_section_vals)
2974 6594 : logger => cp_get_default_logger()
2975 6594 : print_section_vals => section_vals_get_subs_vals(dft_section, "PRINT")
2976 6594 : qs_section_vals => section_vals_get_subs_vals(dft_section, "QS")
2977 6594 : CALL section_vals_get(qs_section_vals, section=qs_section)
2978 :
2979 6594 : NULLIFY (enum, keyword)
2980 6594 : keyword => section_get_keyword(qs_section, "METHOD")
2981 6594 : CALL keyword_get(keyword, enum=enum)
2982 6594 : method = TRIM(enum_i2c(enum, qs_control%method_id))
2983 :
2984 6594 : NULLIFY (enum, keyword)
2985 6594 : keyword => section_get_keyword(qs_section, "QUADRATURE")
2986 6594 : CALL keyword_get(keyword, enum=enum)
2987 6594 : quadrature = TRIM(enum_i2c(enum, qs_control%gapw_control%quadrature))
2988 :
2989 : output_unit = cp_print_key_unit_nr(logger, print_section_vals, &
2990 6594 : "DFT_CONTROL_PARAMETERS", extension=".Log")
2991 6594 : IF (output_unit > 0) THEN
2992 1557 : ngrid_level = SIZE(qs_control%e_cutoff)
2993 : WRITE (UNIT=output_unit, FMT="(/,T2,A,T61,A20)") &
2994 1557 : "QS| Method:", ADJUSTR(method)
2995 1557 : IF (qs_control%pw_grid_opt%spherical) THEN
2996 : WRITE (UNIT=output_unit, FMT="(T2,A,T61,A)") &
2997 0 : "QS| Density plane wave grid type", " SPHERICAL HALFSPACE"
2998 1557 : ELSE IF (qs_control%pw_grid_opt%fullspace) THEN
2999 : WRITE (UNIT=output_unit, FMT="(T2,A,T57,A)") &
3000 1557 : "QS| Density plane wave grid type", " NON-SPHERICAL FULLSPACE"
3001 : ELSE
3002 : WRITE (UNIT=output_unit, FMT="(T2,A,T57,A)") &
3003 0 : "QS| Density plane wave grid type", " NON-SPHERICAL HALFSPACE"
3004 : END IF
3005 : WRITE (UNIT=output_unit, FMT="(T2,A,T71,I10)") &
3006 1557 : "QS| Number of grid levels:", SIZE(qs_control%e_cutoff)
3007 1557 : IF (ngrid_level == 1) THEN
3008 : WRITE (UNIT=output_unit, FMT="(T2,A,T71,F10.1)") &
3009 80 : "QS| Density cutoff [a.u.]:", qs_control%e_cutoff(1)
3010 : ELSE
3011 : WRITE (UNIT=output_unit, FMT="(T2,A,T71,F10.1)") &
3012 1477 : "QS| Density cutoff [a.u.]:", qs_control%cutoff
3013 1477 : IF (qs_control%commensurate_mgrids) THEN
3014 132 : WRITE (UNIT=output_unit, FMT="(T2,A)") "QS| Using commensurate multigrids"
3015 : END IF
3016 : WRITE (UNIT=output_unit, FMT="(T2,A,T71,F10.1)") &
3017 1477 : "QS| Multi grid cutoff [a.u.]: 1) grid level", qs_control%e_cutoff(1)
3018 : WRITE (UNIT=output_unit, FMT="(T2,A,I3,A,T71,F10.1)") &
3019 4600 : ("QS| ", igrid_level, ") grid level", &
3020 6077 : qs_control%e_cutoff(igrid_level), &
3021 7554 : igrid_level=2, SIZE(qs_control%e_cutoff))
3022 : END IF
3023 1557 : IF (qs_control%pao) THEN
3024 0 : WRITE (UNIT=output_unit, FMT="(T2,A)") "QS| PAO active"
3025 : END IF
3026 : WRITE (UNIT=output_unit, FMT="(T2,A,T71,F10.1)") &
3027 1557 : "QS| Grid level progression factor:", qs_control%progression_factor
3028 : WRITE (UNIT=output_unit, FMT="(T2,A,T71,F10.1)") &
3029 1557 : "QS| Relative density cutoff [a.u.]:", qs_control%relative_cutoff
3030 : WRITE (UNIT=output_unit, FMT="(T2,A,T73,ES8.1)") &
3031 1557 : "QS| Interaction thresholds: eps_pgf_orb:", &
3032 1557 : qs_control%eps_pgf_orb, &
3033 1557 : "QS| eps_filter_matrix:", &
3034 1557 : qs_control%eps_filter_matrix, &
3035 1557 : "QS| eps_core_charge:", &
3036 1557 : qs_control%eps_core_charge, &
3037 1557 : "QS| eps_rho_gspace:", &
3038 1557 : qs_control%eps_rho_gspace, &
3039 1557 : "QS| eps_rho_rspace:", &
3040 1557 : qs_control%eps_rho_rspace, &
3041 1557 : "QS| eps_gvg_rspace:", &
3042 1557 : qs_control%eps_gvg_rspace, &
3043 1557 : "QS| eps_ppl:", &
3044 1557 : qs_control%eps_ppl, &
3045 1557 : "QS| eps_ppnl:", &
3046 3114 : qs_control%eps_ppnl
3047 1557 : IF (qs_control%gapw) THEN
3048 293 : IF (qs_control%gapw_control%accurate_xcint) THEN
3049 : WRITE (UNIT=output_unit, FMT="(T2,A,T79,I2)") &
3050 47 : "QS| GAPW| XC integration using accurate scheme: Polynomial order (2n) =", &
3051 94 : qs_control%gapw_control%oweights
3052 : WRITE (UNIT=output_unit, FMT="(T2,A,T69,F12.6)") &
3053 47 : "QS| GAPW| Ref. exponent =", &
3054 94 : qs_control%gapw_control%aweights
3055 : END IF
3056 : !
3057 561 : SELECT CASE (qs_control%gapw_control%basis_1c)
3058 : CASE (gapw_1c_orb)
3059 : WRITE (UNIT=output_unit, FMT="(T2,A)") &
3060 268 : "QS| GAPW| One center basis from orbital basis primitives"
3061 : CASE (gapw_1c_small)
3062 : WRITE (UNIT=output_unit, FMT="(T2,A)") &
3063 21 : "QS| GAPW| One center basis extended with primitives (small:s)"
3064 : CASE (gapw_1c_medium)
3065 : WRITE (UNIT=output_unit, FMT="(T2,A)") &
3066 1 : "QS| GAPW| One center basis extended with primitives (medium:sp)"
3067 : CASE (gapw_1c_large)
3068 : WRITE (UNIT=output_unit, FMT="(T2,A)") &
3069 2 : "QS| GAPW| One center basis extended with primitives (large:spd)"
3070 : CASE (gapw_1c_very_large)
3071 : WRITE (UNIT=output_unit, FMT="(T2,A)") &
3072 1 : "QS| GAPW| One center basis extended with primitives (very large:spdf)"
3073 : CASE DEFAULT
3074 293 : CPABORT("basis_1c incorrect")
3075 : END SELECT
3076 : WRITE (UNIT=output_unit, FMT="(T2,A,T73,ES8.1)") &
3077 293 : "QS| GAPW| eps_fit:", &
3078 293 : qs_control%gapw_control%eps_fit, &
3079 293 : "QS| GAPW| eps_iso:", &
3080 293 : qs_control%gapw_control%eps_iso, &
3081 293 : "QS| GAPW| eps_svd:", &
3082 293 : qs_control%gapw_control%eps_svd, &
3083 293 : "QS| GAPW| eps_cpc:", &
3084 586 : qs_control%gapw_control%eps_cpc
3085 : WRITE (UNIT=output_unit, FMT="(T2,A,T61,A20)") &
3086 293 : "QS| GAPW| atom-r-grid: quadrature:", &
3087 586 : ADJUSTR(quadrature)
3088 : WRITE (UNIT=output_unit, FMT="(T2,A,T71,I10)") &
3089 293 : "QS| GAPW| atom-s-grid: max l :", &
3090 293 : qs_control%gapw_control%lmax_sphere, &
3091 293 : "QS| GAPW| max_l_rho0 :", &
3092 586 : qs_control%gapw_control%lmax_rho0
3093 293 : IF (qs_control%gapw_control%non_paw_atoms) THEN
3094 : WRITE (UNIT=output_unit, FMT="(T2,A)") &
3095 62 : "QS| GAPW| At least one kind is NOT PAW, i.e. it has only soft AO "
3096 : END IF
3097 293 : IF (qs_control%gapw_control%nopaw_as_gpw) THEN
3098 : WRITE (UNIT=output_unit, FMT="(T2,A)") &
3099 62 : "QS| GAPW| The NOT PAW atoms are treated fully GPW"
3100 : END IF
3101 : END IF
3102 1557 : IF (qs_control%gapw_xc) THEN
3103 71 : SELECT CASE (qs_control%gapw_control%basis_1c)
3104 : CASE (gapw_1c_orb)
3105 : WRITE (UNIT=output_unit, FMT="(T2,A)") &
3106 32 : "QS| GAPW_XC| One center basis from orbital basis primitives"
3107 : CASE (gapw_1c_small)
3108 : WRITE (UNIT=output_unit, FMT="(T2,A)") &
3109 7 : "QS| GAPW_XC| One center basis extended with primitives (small:s)"
3110 : CASE (gapw_1c_medium)
3111 : WRITE (UNIT=output_unit, FMT="(T2,A)") &
3112 0 : "QS| GAPW_XC| One center basis extended with primitives (medium:sp)"
3113 : CASE (gapw_1c_large)
3114 : WRITE (UNIT=output_unit, FMT="(T2,A)") &
3115 0 : "QS| GAPW_XC| One center basis extended with primitives (large:spd)"
3116 : CASE (gapw_1c_very_large)
3117 : WRITE (UNIT=output_unit, FMT="(T2,A)") &
3118 0 : "QS| GAPW_XC| One center basis extended with primitives (very large:spdf)"
3119 : CASE DEFAULT
3120 39 : CPABORT("basis_1c incorrect")
3121 : END SELECT
3122 : WRITE (UNIT=output_unit, FMT="(T2,A,T73,ES8.1)") &
3123 39 : "QS| GAPW_XC| eps_fit:", &
3124 39 : qs_control%gapw_control%eps_fit, &
3125 39 : "QS| GAPW_XC| eps_iso:", &
3126 39 : qs_control%gapw_control%eps_iso, &
3127 39 : "QS| GAPW_XC| eps_svd:", &
3128 78 : qs_control%gapw_control%eps_svd
3129 : WRITE (UNIT=output_unit, FMT="(T2,A,T55,A30)") &
3130 39 : "QS| GAPW_XC|atom-r-grid: quadrature:", &
3131 78 : enum_i2c(enum, qs_control%gapw_control%quadrature)
3132 : WRITE (UNIT=output_unit, FMT="(T2,A,T71,I10)") &
3133 39 : "QS| GAPW_XC| atom-s-grid: max l :", &
3134 78 : qs_control%gapw_control%lmax_sphere
3135 : END IF
3136 1557 : IF (qs_control%mulliken_restraint) THEN
3137 : WRITE (UNIT=output_unit, FMT="(T2,A,T73,ES8.1)") &
3138 1 : "QS| Mulliken restraint target", qs_control%mulliken_restraint_control%target
3139 : WRITE (UNIT=output_unit, FMT="(T2,A,T73,ES8.1)") &
3140 1 : "QS| Mulliken restraint strength", qs_control%mulliken_restraint_control%strength
3141 : WRITE (UNIT=output_unit, FMT="(T2,A,T73,I8)") &
3142 1 : "QS| Mulliken restraint atoms: ", qs_control%mulliken_restraint_control%natoms
3143 2 : WRITE (UNIT=output_unit, FMT="(5I8)") qs_control%mulliken_restraint_control%atoms
3144 : END IF
3145 1557 : IF (qs_control%ddapc_restraint) THEN
3146 14 : DO i = 1, SIZE(qs_control%ddapc_restraint_control)
3147 8 : ddapc_restraint_control => qs_control%ddapc_restraint_control(i)
3148 8 : IF (SIZE(qs_control%ddapc_restraint_control) > 1) THEN
3149 : WRITE (UNIT=output_unit, FMT="(T2,A,T3,I8)") &
3150 3 : "QS| parameters for DDAPC restraint number", i
3151 : END IF
3152 : WRITE (UNIT=output_unit, FMT="(T2,A,T73,ES8.1)") &
3153 8 : "QS| ddapc restraint target", ddapc_restraint_control%target
3154 : WRITE (UNIT=output_unit, FMT="(T2,A,T73,ES8.1)") &
3155 8 : "QS| ddapc restraint strength", ddapc_restraint_control%strength
3156 : WRITE (UNIT=output_unit, FMT="(T2,A,T73,I8)") &
3157 8 : "QS| ddapc restraint atoms: ", ddapc_restraint_control%natoms
3158 17 : WRITE (UNIT=output_unit, FMT="(5I8)") ddapc_restraint_control%atoms
3159 8 : WRITE (UNIT=output_unit, FMT="(T2,A)") "Coefficients:"
3160 17 : WRITE (UNIT=output_unit, FMT="(5F6.2)") ddapc_restraint_control%coeff
3161 6 : SELECT CASE (ddapc_restraint_control%functional_form)
3162 : CASE (do_ddapc_restraint)
3163 : WRITE (UNIT=output_unit, FMT="(T2,A,T61,A20)") &
3164 3 : "QS| ddapc restraint functional form :", "RESTRAINT"
3165 : CASE (do_ddapc_constraint)
3166 : WRITE (UNIT=output_unit, FMT="(T2,A,T61,A20)") &
3167 5 : "QS| ddapc restraint functional form :", "CONSTRAINT"
3168 : CASE DEFAULT
3169 8 : CPABORT("Unknown ddapc restraint")
3170 : END SELECT
3171 : END DO
3172 : END IF
3173 1557 : IF (qs_control%s2_restraint) THEN
3174 : WRITE (UNIT=output_unit, FMT="(T2,A,T73,ES8.1)") &
3175 0 : "QS| s2 restraint target", qs_control%s2_restraint_control%target
3176 : WRITE (UNIT=output_unit, FMT="(T2,A,T73,ES8.1)") &
3177 0 : "QS| s2 restraint strength", qs_control%s2_restraint_control%strength
3178 0 : SELECT CASE (qs_control%s2_restraint_control%functional_form)
3179 : CASE (do_s2_restraint)
3180 : WRITE (UNIT=output_unit, FMT="(T2,A,T61,A20)") &
3181 0 : "QS| s2 restraint functional form :", "RESTRAINT"
3182 0 : CPABORT("Not yet implemented")
3183 : CASE (do_s2_constraint)
3184 : WRITE (UNIT=output_unit, FMT="(T2,A,T61,A20)") &
3185 0 : "QS| s2 restraint functional form :", "CONSTRAINT"
3186 : CASE DEFAULT
3187 0 : CPABORT("Unknown ddapc restraint")
3188 : END SELECT
3189 : END IF
3190 : END IF
3191 : CALL cp_print_key_finished_output(output_unit, logger, print_section_vals, &
3192 6594 : "DFT_CONTROL_PARAMETERS")
3193 :
3194 6594 : CALL timestop(handle)
3195 :
3196 : END SUBROUTINE write_qs_control
3197 :
3198 : ! **************************************************************************************************
3199 : !> \brief reads the input parameters needed for ddapc.
3200 : !> \param qs_control ...
3201 : !> \param qs_section ...
3202 : !> \param ddapc_restraint_section ...
3203 : !> \author fschiff
3204 : !> \note
3205 : !> either reads DFT%QS%DDAPC_RESTRAINT or PROPERTIES%ET_coupling
3206 : !> if(qs_section is present the DFT part is read, if ddapc_restraint_section
3207 : !> is present ET_COUPLING is read. Avoid having both!!!
3208 : ! **************************************************************************************************
3209 14 : SUBROUTINE read_ddapc_section(qs_control, qs_section, ddapc_restraint_section)
3210 :
3211 : TYPE(qs_control_type), INTENT(INOUT) :: qs_control
3212 : TYPE(section_vals_type), OPTIONAL, POINTER :: qs_section, ddapc_restraint_section
3213 :
3214 : INTEGER :: i, j, jj, k, n_rep
3215 14 : INTEGER, DIMENSION(:), POINTER :: tmplist
3216 14 : REAL(KIND=dp), DIMENSION(:), POINTER :: rtmplist
3217 : TYPE(ddapc_restraint_type), POINTER :: ddapc_restraint_control
3218 : TYPE(section_vals_type), POINTER :: ddapc_section
3219 :
3220 14 : IF (PRESENT(ddapc_restraint_section)) THEN
3221 0 : IF (ASSOCIATED(qs_control%ddapc_restraint_control)) THEN
3222 0 : IF (SIZE(qs_control%ddapc_restraint_control) >= 2) THEN
3223 0 : CPABORT("ET_COUPLING cannot be used in combination with a normal restraint")
3224 : END IF
3225 : ELSE
3226 0 : ddapc_section => ddapc_restraint_section
3227 0 : ALLOCATE (qs_control%ddapc_restraint_control(1))
3228 : END IF
3229 : END IF
3230 :
3231 14 : IF (PRESENT(qs_section)) THEN
3232 14 : NULLIFY (ddapc_section)
3233 : ddapc_section => section_vals_get_subs_vals(qs_section, &
3234 14 : "DDAPC_RESTRAINT")
3235 : END IF
3236 :
3237 32 : DO i = 1, SIZE(qs_control%ddapc_restraint_control)
3238 :
3239 18 : CALL ddapc_control_create(qs_control%ddapc_restraint_control(i))
3240 18 : ddapc_restraint_control => qs_control%ddapc_restraint_control(i)
3241 :
3242 : CALL section_vals_val_get(ddapc_section, "STRENGTH", i_rep_section=i, &
3243 18 : r_val=ddapc_restraint_control%strength)
3244 : CALL section_vals_val_get(ddapc_section, "TARGET", i_rep_section=i, &
3245 18 : r_val=ddapc_restraint_control%target)
3246 : CALL section_vals_val_get(ddapc_section, "FUNCTIONAL_FORM", i_rep_section=i, &
3247 18 : i_val=ddapc_restraint_control%functional_form)
3248 : CALL section_vals_val_get(ddapc_section, "ATOMS", i_rep_section=i, &
3249 18 : n_rep_val=n_rep)
3250 : CALL section_vals_val_get(ddapc_section, "TYPE_OF_DENSITY", i_rep_section=i, &
3251 18 : i_val=ddapc_restraint_control%density_type)
3252 :
3253 18 : jj = 0
3254 36 : DO k = 1, n_rep
3255 : CALL section_vals_val_get(ddapc_section, "ATOMS", i_rep_section=i, &
3256 18 : i_rep_val=k, i_vals=tmplist)
3257 56 : DO j = 1, SIZE(tmplist)
3258 38 : jj = jj + 1
3259 : END DO
3260 : END DO
3261 18 : IF (jj < 1) CPABORT("Need at least 1 atom to use ddapc constraints")
3262 18 : ddapc_restraint_control%natoms = jj
3263 18 : IF (ASSOCIATED(ddapc_restraint_control%atoms)) THEN
3264 0 : DEALLOCATE (ddapc_restraint_control%atoms)
3265 : END IF
3266 54 : ALLOCATE (ddapc_restraint_control%atoms(ddapc_restraint_control%natoms))
3267 18 : jj = 0
3268 36 : DO k = 1, n_rep
3269 : CALL section_vals_val_get(ddapc_section, "ATOMS", i_rep_section=i, &
3270 18 : i_rep_val=k, i_vals=tmplist)
3271 56 : DO j = 1, SIZE(tmplist)
3272 20 : jj = jj + 1
3273 38 : ddapc_restraint_control%atoms(jj) = tmplist(j)
3274 : END DO
3275 : END DO
3276 :
3277 18 : IF (ASSOCIATED(ddapc_restraint_control%coeff)) THEN
3278 0 : DEALLOCATE (ddapc_restraint_control%coeff)
3279 : END IF
3280 54 : ALLOCATE (ddapc_restraint_control%coeff(ddapc_restraint_control%natoms))
3281 38 : ddapc_restraint_control%coeff = 1.0_dp
3282 :
3283 : CALL section_vals_val_get(ddapc_section, "COEFF", i_rep_section=i, &
3284 18 : n_rep_val=n_rep)
3285 18 : jj = 0
3286 20 : DO k = 1, n_rep
3287 : CALL section_vals_val_get(ddapc_section, "COEFF", i_rep_section=i, &
3288 2 : i_rep_val=k, r_vals=rtmplist)
3289 22 : DO j = 1, SIZE(rtmplist)
3290 2 : jj = jj + 1
3291 2 : IF (jj > ddapc_restraint_control%natoms) THEN
3292 0 : CPABORT("Need the same number of coeff as there are atoms ")
3293 : END IF
3294 4 : ddapc_restraint_control%coeff(jj) = rtmplist(j)
3295 : END DO
3296 : END DO
3297 68 : IF (jj < ddapc_restraint_control%natoms .AND. jj /= 0) THEN
3298 0 : CPABORT("Need no or the same number of coeff as there are atoms.")
3299 : END IF
3300 : END DO
3301 14 : k = 0
3302 32 : DO i = 1, SIZE(qs_control%ddapc_restraint_control)
3303 18 : IF (qs_control%ddapc_restraint_control(i)%functional_form == &
3304 24 : do_ddapc_constraint) k = k + 1
3305 : END DO
3306 14 : IF (k == 2) CALL cp_abort(__LOCATION__, &
3307 0 : "Only a single constraint possible yet, try to use restraints instead ")
3308 :
3309 14 : END SUBROUTINE read_ddapc_section
3310 :
3311 : ! **************************************************************************************************
3312 : !> \brief ...
3313 : !> \param dft_control ...
3314 : !> \param efield_section ...
3315 : !> \param cell ...
3316 : ! **************************************************************************************************
3317 342 : SUBROUTINE read_efield_sections(dft_control, efield_section, cell)
3318 : TYPE(dft_control_type), POINTER :: dft_control
3319 : TYPE(section_vals_type), POINTER :: efield_section
3320 : TYPE(cell_type), OPTIONAL, POINTER :: cell
3321 :
3322 : CHARACTER(len=default_path_length) :: file_name
3323 : INTEGER :: i, io, j, n, unit_nr
3324 : LOGICAL :: amplitude_explicit, intensity_explicit
3325 342 : REAL(KIND=dp), DIMENSION(:), POINTER :: tmp_vals
3326 : TYPE(efield_type), POINTER :: efield
3327 : TYPE(section_vals_type), POINTER :: tmp_section
3328 :
3329 684 : DO i = 1, SIZE(dft_control%efield_fields)
3330 342 : NULLIFY (dft_control%efield_fields(i)%efield)
3331 1368 : ALLOCATE (dft_control%efield_fields(i)%efield)
3332 342 : efield => dft_control%efield_fields(i)%efield
3333 342 : NULLIFY (efield%envelop_i_vars, efield%envelop_r_vars)
3334 : CALL section_vals_val_get(efield_section, "INTENSITY", i_rep_section=i, &
3335 342 : r_val=efield%strength, explicit=intensity_explicit)
3336 : CALL section_vals_val_get(efield_section, "AMPLITUDE", i_rep_section=i, &
3337 342 : r_val=efield%amplitude, explicit=amplitude_explicit)
3338 :
3339 342 : IF (intensity_explicit .AND. amplitude_explicit) THEN
3340 0 : CPABORT("Both INTENSITY and AMPLITUDE provided in EFIELD section.")
3341 : END IF
3342 :
3343 342 : IF (intensity_explicit) THEN
3344 48 : efield%amplitude = SQRT(efield%strength/(3.50944_dp*10.0_dp**16))
3345 : END IF
3346 :
3347 342 : IF (amplitude_explicit) THEN
3348 0 : efield%strength = (3.50944_dp*10.0_dp**16)*(efield%amplitude)**2
3349 : END IF
3350 :
3351 : CALL section_vals_val_get(efield_section, "POLARISATION", i_rep_section=i, &
3352 342 : r_vals=tmp_vals)
3353 1026 : ALLOCATE (efield%polarisation(SIZE(tmp_vals)))
3354 2394 : efield%polarisation = tmp_vals
3355 342 : IF (PRESENT(cell)) THEN
3356 342 : IF (ASSOCIATED(cell)) CALL cell_transform_input_cartesian(cell, efield%polarisation(1:3))
3357 : END IF
3358 : CALL section_vals_val_get(efield_section, "PHASE", i_rep_section=i, &
3359 342 : r_val=efield%phase_offset)
3360 : CALL section_vals_val_get(efield_section, "ENVELOP", i_rep_section=i, &
3361 342 : i_val=efield%envelop_id)
3362 : CALL section_vals_val_get(efield_section, "WAVELENGTH", i_rep_section=i, &
3363 342 : r_val=efield%wavelength)
3364 : CALL section_vals_val_get(efield_section, "VEC_POT_INITIAL", i_rep_section=i, &
3365 342 : r_vals=tmp_vals)
3366 2394 : efield%vec_pot_initial = tmp_vals
3367 342 : IF (PRESENT(cell)) THEN
3368 342 : IF (ASSOCIATED(cell)) CALL cell_transform_input_cartesian(cell, efield%vec_pot_initial(1:3))
3369 : END IF
3370 :
3371 1026 : IF (efield%envelop_id == constant_env) THEN
3372 326 : ALLOCATE (efield%envelop_i_vars(2))
3373 326 : tmp_section => section_vals_get_subs_vals(efield_section, "CONSTANT_ENV", i_rep_section=i)
3374 : CALL section_vals_val_get(tmp_section, "START_STEP", &
3375 326 : i_val=efield%envelop_i_vars(1))
3376 : CALL section_vals_val_get(tmp_section, "END_STEP", &
3377 326 : i_val=efield%envelop_i_vars(2))
3378 16 : ELSE IF (efield%envelop_id == gaussian_env) THEN
3379 12 : ALLOCATE (efield%envelop_r_vars(2))
3380 12 : tmp_section => section_vals_get_subs_vals(efield_section, "GAUSSIAN_ENV", i_rep_section=i)
3381 : CALL section_vals_val_get(tmp_section, "T0", &
3382 12 : r_val=efield%envelop_r_vars(1))
3383 : CALL section_vals_val_get(tmp_section, "SIGMA", &
3384 12 : r_val=efield%envelop_r_vars(2))
3385 4 : ELSE IF (efield%envelop_id == ramp_env) THEN
3386 2 : ALLOCATE (efield%envelop_i_vars(4))
3387 2 : tmp_section => section_vals_get_subs_vals(efield_section, "RAMP_ENV", i_rep_section=i)
3388 : CALL section_vals_val_get(tmp_section, "START_STEP_IN", &
3389 2 : i_val=efield%envelop_i_vars(1))
3390 : CALL section_vals_val_get(tmp_section, "END_STEP_IN", &
3391 2 : i_val=efield%envelop_i_vars(2))
3392 : CALL section_vals_val_get(tmp_section, "START_STEP_OUT", &
3393 2 : i_val=efield%envelop_i_vars(3))
3394 : CALL section_vals_val_get(tmp_section, "END_STEP_OUT", &
3395 2 : i_val=efield%envelop_i_vars(4))
3396 2 : ELSE IF (efield%envelop_id == custom_env) THEN
3397 2 : tmp_section => section_vals_get_subs_vals(efield_section, "CUSTOM_ENV", i_rep_section=i)
3398 2 : CALL section_vals_val_get(tmp_section, "EFIELD_FILE_NAME", c_val=file_name)
3399 2 : CALL open_file(file_name=TRIM(file_name), file_action="READ", file_status="OLD", unit_number=unit_nr)
3400 : !Determine the number of lines in file
3401 2 : n = 0
3402 10 : DO WHILE (.TRUE.)
3403 12 : READ (unit_nr, *, iostat=io)
3404 12 : IF (io /= 0) EXIT
3405 10 : n = n + 1
3406 : END DO
3407 2 : REWIND (unit_nr)
3408 6 : ALLOCATE (efield%envelop_r_vars(n + 1))
3409 : !Store the timestep of the list in the first entry of the r_vars
3410 2 : CALL section_vals_val_get(tmp_section, "TIMESTEP", r_val=efield%envelop_r_vars(1))
3411 : !Read the file
3412 12 : DO j = 2, n + 1
3413 10 : READ (unit_nr, *) efield%envelop_r_vars(j)
3414 12 : efield%envelop_r_vars(j) = cp_unit_to_cp2k(efield%envelop_r_vars(j), "volt/m")
3415 : END DO
3416 2 : CALL close_file(unit_nr)
3417 : END IF
3418 : END DO
3419 342 : END SUBROUTINE read_efield_sections
3420 :
3421 : ! **************************************************************************************************
3422 : !> \brief reads the input parameters needed real time propagation
3423 : !> \param dft_control ...
3424 : !> \param rtp_section ...
3425 : !> \author fschiff
3426 : ! **************************************************************************************************
3427 2268 : SUBROUTINE read_rtp_section(dft_control, rtp_section)
3428 :
3429 : TYPE(dft_control_type), INTENT(INOUT) :: dft_control
3430 : TYPE(section_vals_type), POINTER :: rtp_section
3431 :
3432 : INTEGER :: i, j, n_elems
3433 324 : INTEGER, DIMENSION(:), POINTER :: tmp
3434 : LOGICAL :: is_present, linearize_bse_propagation, &
3435 : local_moment_possible
3436 : TYPE(section_vals_type), POINTER :: proj_mo_section, subsection
3437 :
3438 3888 : ALLOCATE (dft_control%rtp_control)
3439 : CALL section_vals_val_get(rtp_section, "MAX_ITER", &
3440 324 : i_val=dft_control%rtp_control%max_iter)
3441 : CALL section_vals_val_get(rtp_section, "MAT_EXP", &
3442 324 : i_val=dft_control%rtp_control%mat_exp)
3443 : CALL section_vals_val_get(rtp_section, "ASPC_ORDER", &
3444 324 : i_val=dft_control%rtp_control%aspc_order)
3445 : CALL section_vals_val_get(rtp_section, "EXP_ACCURACY", &
3446 324 : r_val=dft_control%rtp_control%eps_exp)
3447 : CALL section_vals_val_get(rtp_section, "RTBSE%_SECTION_PARAMETERS_", &
3448 324 : i_val=dft_control%rtp_control%rtp_method)
3449 : CALL section_vals_val_get(rtp_section, "RTBSE%RTBSE_HAMILTONIAN", &
3450 324 : i_val=dft_control%rtp_control%rtbse_ham)
3451 : CALL section_vals_val_get(rtp_section, "RTBSE%LINEARIZED_BSE_PROPAGATION", &
3452 324 : l_val=linearize_bse_propagation)
3453 : ! Change rtp_method to linearized bse. The section parameter also feeds bs_env%rtp_method,
3454 : ! which gates the W(w=0) build in the GW step - TDDFT there would dispatch the linearized
3455 : ! propagator with no screened interaction to propagate with, so reject the combination.
3456 324 : IF (linearize_bse_propagation) THEN
3457 58 : IF (dft_control%rtp_control%rtp_method /= rtp_method_bse) THEN
3458 : CALL cp_abort(__LOCATION__, &
3459 : "LINEARIZED_BSE_PROPAGATION requires the RTBSE section opened as "// &
3460 0 : "'&RTBSE' or '&RTBSE RTBSE', not '&RTBSE TDDFT'.")
3461 : END IF
3462 58 : dft_control%rtp_control%rtp_method = rtp_method_bse_linearized
3463 : END IF
3464 :
3465 : CALL section_vals_val_get(rtp_section, "PROPAGATOR", &
3466 324 : i_val=dft_control%rtp_control%propagator)
3467 : CALL section_vals_val_get(rtp_section, "EPS_ITER", &
3468 324 : r_val=dft_control%rtp_control%eps_ener)
3469 : CALL section_vals_val_get(rtp_section, "INITIAL_WFN", &
3470 324 : i_val=dft_control%rtp_control%initial_wfn)
3471 : CALL section_vals_val_get(rtp_section, "HFX_BALANCE_IN_CORE", &
3472 324 : l_val=dft_control%rtp_control%hfx_redistribute)
3473 : CALL section_vals_val_get(rtp_section, "APPLY_WFN_MIX_INIT_RESTART", &
3474 324 : l_val=dft_control%rtp_control%apply_wfn_mix_init_restart)
3475 : CALL section_vals_val_get(rtp_section, "APPLY_DELTA_PULSE", &
3476 324 : l_val=dft_control%rtp_control%apply_delta_pulse)
3477 : CALL section_vals_val_get(rtp_section, "APPLY_DELTA_PULSE_MAG", &
3478 324 : l_val=dft_control%rtp_control%apply_delta_pulse_mag)
3479 : CALL section_vals_val_get(rtp_section, "VELOCITY_GAUGE", &
3480 324 : l_val=dft_control%rtp_control%velocity_gauge)
3481 : CALL section_vals_val_get(rtp_section, "VG_COM_NL", &
3482 324 : l_val=dft_control%rtp_control%nl_gauge_transform)
3483 : CALL section_vals_val_get(rtp_section, "PERIODIC", &
3484 324 : l_val=dft_control%rtp_control%periodic)
3485 : CALL section_vals_val_get(rtp_section, "DENSITY_PROPAGATION", &
3486 324 : l_val=dft_control%rtp_control%linear_scaling)
3487 : CALL section_vals_val_get(rtp_section, "MCWEENY_MAX_ITER", &
3488 324 : i_val=dft_control%rtp_control%mcweeny_max_iter)
3489 : CALL section_vals_val_get(rtp_section, "ACCURACY_REFINEMENT", &
3490 324 : i_val=dft_control%rtp_control%acc_ref)
3491 : CALL section_vals_val_get(rtp_section, "MCWEENY_EPS", &
3492 324 : r_val=dft_control%rtp_control%mcweeny_eps)
3493 : CALL section_vals_val_get(rtp_section, "DELTA_PULSE_SCALE", &
3494 324 : r_val=dft_control%rtp_control%delta_pulse_scale)
3495 : CALL section_vals_val_get(rtp_section, "DELTA_PULSE_DIRECTION", &
3496 324 : i_vals=tmp)
3497 1296 : dft_control%rtp_control%delta_pulse_direction = tmp
3498 : CALL section_vals_val_get(rtp_section, "SC_CHECK_START", &
3499 324 : i_val=dft_control%rtp_control%sc_check_start)
3500 324 : proj_mo_section => section_vals_get_subs_vals(rtp_section, "PRINT%PROJECTION_MO")
3501 324 : CALL section_vals_get(proj_mo_section, explicit=is_present)
3502 324 : IF (is_present) THEN
3503 4 : IF (dft_control%rtp_control%linear_scaling) THEN
3504 : CALL cp_abort(__LOCATION__, &
3505 : "You have defined a time dependent projection of mos, but "// &
3506 : "only the density matrix is propagated (DENSITY_PROPAGATION "// &
3507 : ".TRUE.). Please either use MO-based real time DFT or do not "// &
3508 0 : "define any PRINT%PROJECTION_MO section")
3509 : END IF
3510 4 : dft_control%rtp_control%is_proj_mo = .TRUE.
3511 : ELSE
3512 320 : dft_control%rtp_control%is_proj_mo = .FALSE.
3513 : END IF
3514 : ! Moment trace
3515 : local_moment_possible = (dft_control%rtp_control%rtp_method == rtp_method_bse .OR. &
3516 : dft_control%rtp_control%rtp_method == rtp_method_bse_linearized) .OR. &
3517 324 : ((.NOT. dft_control%rtp_control%periodic) .AND. dft_control%rtp_control%linear_scaling)
3518 : ! TODO : Implement for other moment operators
3519 324 : subsection => section_vals_get_subs_vals(rtp_section, "PRINT%MOMENTS")
3520 324 : CALL section_vals_get(subsection, explicit=is_present)
3521 : ! Trigger the flag
3522 : dft_control%rtp_control%save_local_moments = &
3523 324 : is_present .OR. dft_control%rtp_control%save_local_moments
3524 324 : IF (is_present .AND. (.NOT. local_moment_possible)) THEN
3525 : CALL cp_abort(__LOCATION__, "Moments trace printing only "// &
3526 : "implemented in non-periodic systems in linear scaling. "// &
3527 0 : "Please use DFT%PRINT%MOMENTS for other printing.")
3528 : END IF
3529 : CALL section_vals_val_get(rtp_section, "PRINT%MOMENTS%REFERENCE", &
3530 324 : i_val=dft_control%rtp_control%moment_trace_ref_type)
3531 : CALL section_vals_val_get(rtp_section, "PRINT%MOMENTS%REFERENCE_POINT", &
3532 324 : r_vals=dft_control%rtp_control%moment_trace_user_ref_point)
3533 : ! Moment Fourier transform
3534 324 : subsection => section_vals_get_subs_vals(rtp_section, "PRINT%MOMENTS_FT")
3535 324 : CALL section_vals_get(subsection, explicit=is_present)
3536 : ! Trigger the flag
3537 : dft_control%rtp_control%save_local_moments = &
3538 324 : is_present .OR. dft_control%rtp_control%save_local_moments
3539 324 : IF (is_present .AND. (.NOT. local_moment_possible)) THEN
3540 : ! Not implemented
3541 : CALL cp_abort(__LOCATION__, "Moments Fourier transform printing "// &
3542 0 : "implemented only for non-periodic systems in linear scaling.")
3543 : END IF
3544 : ! General FT settings
3545 : CALL section_vals_val_get(rtp_section, "FT%DAMPING", &
3546 324 : r_val=dft_control%rtp_control%ft_damping)
3547 : CALL section_vals_val_get(rtp_section, "FT%START_TIME", &
3548 324 : r_val=dft_control%rtp_control%ft_t0)
3549 : ! Padé settings
3550 324 : subsection => section_vals_get_subs_vals(rtp_section, "FT%PADE")
3551 : CALL section_vals_val_get(subsection, "_SECTION_PARAMETERS_", &
3552 324 : l_val=dft_control%rtp_control%pade_requested)
3553 : CALL section_vals_val_get(subsection, "E_MIN", &
3554 324 : r_val=dft_control%rtp_control%pade_e_min)
3555 : CALL section_vals_val_get(subsection, "E_STEP", &
3556 324 : r_val=dft_control%rtp_control%pade_e_step)
3557 : CALL section_vals_val_get(subsection, "E_MAX", &
3558 324 : r_val=dft_control%rtp_control%pade_e_max)
3559 : CALL section_vals_val_get(subsection, "FIT_E_MIN", &
3560 324 : r_val=dft_control%rtp_control%pade_fit_e_min)
3561 : CALL section_vals_val_get(subsection, "FIT_E_MAX", &
3562 324 : r_val=dft_control%rtp_control%pade_fit_e_max)
3563 : ! If default settings used for fit_e_min/max, rewrite with appropriate values
3564 324 : IF (dft_control%rtp_control%pade_fit_e_min < 0) THEN
3565 324 : dft_control%rtp_control%pade_fit_e_min = dft_control%rtp_control%pade_e_min
3566 : END IF
3567 324 : IF (dft_control%rtp_control%pade_fit_e_max < 0) THEN
3568 324 : dft_control%rtp_control%pade_fit_e_max = dft_control%rtp_control%pade_e_max
3569 : END IF
3570 : ! Polarizability settings
3571 324 : subsection => section_vals_get_subs_vals(rtp_section, "PRINT%POLARIZABILITY")
3572 324 : CALL section_vals_get(subsection, explicit=is_present)
3573 : ! Trigger the flag
3574 : dft_control%rtp_control%save_local_moments = &
3575 324 : is_present .OR. dft_control%rtp_control%save_local_moments
3576 324 : IF (is_present .AND. (.NOT. local_moment_possible)) THEN
3577 : ! Not implemented
3578 : CALL cp_abort(__LOCATION__, "Polarizability printing "// &
3579 0 : "implemented only for non-periodic systems.")
3580 : END IF
3581 324 : CALL section_vals_val_get(subsection, "ELEMENT", explicit=is_present, n_rep_val=n_elems)
3582 324 : NULLIFY (dft_control%rtp_control%print_pol_elements)
3583 324 : IF (is_present) THEN
3584 : ! Explicit list of elements
3585 : ! Allocate the array
3586 0 : ALLOCATE (dft_control%rtp_control%print_pol_elements(n_elems, 2))
3587 0 : DO i = 1, n_elems
3588 0 : CALL section_vals_val_get(subsection, "ELEMENT", i_vals=tmp, i_rep_val=i)
3589 0 : dft_control%rtp_control%print_pol_elements(i, :) = tmp(:)
3590 : END DO
3591 : ! Do basic sanity checks for pol_element
3592 0 : DO i = 1, n_elems
3593 0 : DO j = 1, 2
3594 0 : IF (dft_control%rtp_control%print_pol_elements(i, j) > 3 .OR. &
3595 0 : dft_control%rtp_control%print_pol_elements(i, j) < 1) THEN
3596 0 : CPABORT("Polarisation tensor element not 1,2 or 3 in at least one index")
3597 : END IF
3598 : END DO
3599 : END DO
3600 : END IF
3601 :
3602 : ! Finally, allow printing of FT observables also in the case when they are not explicitly
3603 : ! required, but they are available, i.e. non-periodic linear scaling calculation
3604 : dft_control%rtp_control%save_local_moments = &
3605 : dft_control%rtp_control%save_local_moments .OR. &
3606 324 : ((.NOT. dft_control%rtp_control%periodic) .AND. dft_control%rtp_control%linear_scaling)
3607 :
3608 324 : END SUBROUTINE read_rtp_section
3609 : ! **************************************************************************************************
3610 : !> \brief Tries to guess the elements of polarization to print
3611 : !> \param dftc DFT parameters
3612 : !> \param elems 2D array, where the guessed element indeces are stored
3613 : !> \date 11.2025
3614 : !> \author Stepan Marek
3615 : ! **************************************************************************************************
3616 90 : SUBROUTINE guess_pol_elements(dftc, elems)
3617 : TYPE(dft_control_type) :: dftc
3618 : INTEGER, DIMENSION(:, :), POINTER :: elems
3619 :
3620 : INTEGER :: i, i_nonzero, n_nonzero
3621 : LOGICAL :: pol_vector_known
3622 : REAL(kind=dp), DIMENSION(3) :: pol_vector
3623 :
3624 90 : pol_vector_known = .FALSE.
3625 :
3626 : ! TODO : More relevant elements for magnetic pulse?
3627 90 : IF (dftc%rtp_control%apply_delta_pulse .OR. dftc%rtp_control%apply_delta_pulse_mag) THEN
3628 336 : pol_vector(:) = REAL(dftc%rtp_control%delta_pulse_direction(:), kind=dp)
3629 : ELSE
3630 : ! Maybe RT field is applied?
3631 24 : pol_vector(:) = dftc%efield_fields(1)%efield%polarisation(:)
3632 : END IF
3633 360 : IF (DOT_PRODUCT(pol_vector, pol_vector) > 0.0_dp) pol_vector_known = .TRUE.
3634 :
3635 : IF (.NOT. pol_vector_known) THEN
3636 0 : CPABORT("Cannot guess polarization elements - please specify!")
3637 : ELSE
3638 : ! Check whether just one element is non-zero
3639 : n_nonzero = 0
3640 360 : DO i = 1, 3
3641 360 : IF (pol_vector(i) /= 0.0_dp) THEN
3642 90 : n_nonzero = n_nonzero + 1
3643 90 : i_nonzero = i
3644 : END IF
3645 : END DO
3646 90 : IF (n_nonzero > 1) THEN
3647 : CALL cp_abort(__LOCATION__, &
3648 : "More than one non-zero field elements - "// &
3649 0 : "cannot guess polarizability elements - please specify!")
3650 90 : ELSE IF (n_nonzero == 0) THEN
3651 : CALL cp_abort(__LOCATION__, &
3652 : "No non-zero field elements - "// &
3653 0 : "cannot guess polarizability elements - please specify!")
3654 : ELSE
3655 : ! Clear guess can be made
3656 : NULLIFY (elems)
3657 90 : ALLOCATE (elems(3, 2))
3658 360 : DO i = 1, 3
3659 270 : elems(i, 1) = i
3660 360 : elems(i, 2) = i_nonzero
3661 : END DO
3662 : END IF
3663 : END IF
3664 90 : END SUBROUTINE guess_pol_elements
3665 :
3666 : ! **************************************************************************************************
3667 : !> \brief Parses the BLOCK_LIST keywords from the ADMM section
3668 : !> \param admm_control ...
3669 : !> \param dft_section ...
3670 : ! **************************************************************************************************
3671 524 : SUBROUTINE read_admm_block_list(admm_control, dft_section)
3672 : TYPE(admm_control_type), POINTER :: admm_control
3673 : TYPE(section_vals_type), POINTER :: dft_section
3674 :
3675 : INTEGER :: irep, list_size, n_rep
3676 524 : INTEGER, DIMENSION(:), POINTER :: tmplist
3677 :
3678 524 : NULLIFY (tmplist)
3679 :
3680 : CALL section_vals_val_get(dft_section, "AUXILIARY_DENSITY_MATRIX_METHOD%BLOCK_LIST", &
3681 524 : n_rep_val=n_rep)
3682 :
3683 1102 : ALLOCATE (admm_control%blocks(n_rep))
3684 :
3685 560 : DO irep = 1, n_rep
3686 : CALL section_vals_val_get(dft_section, "AUXILIARY_DENSITY_MATRIX_METHOD%BLOCK_LIST", &
3687 36 : i_rep_val=irep, i_vals=tmplist)
3688 36 : list_size = SIZE(tmplist)
3689 108 : ALLOCATE (admm_control%blocks(irep)%list(list_size))
3690 732 : admm_control%blocks(irep)%list(:) = tmplist(:)
3691 : END DO
3692 :
3693 524 : END SUBROUTINE read_admm_block_list
3694 :
3695 : ! **************************************************************************************************
3696 : !> \brief ...
3697 : !> \param dft_control ...
3698 : !> \param hairy_probes_section ...
3699 : !> \param
3700 : !> \param
3701 : ! **************************************************************************************************
3702 4 : SUBROUTINE read_hairy_probes_sections(dft_control, hairy_probes_section)
3703 : TYPE(dft_control_type), POINTER :: dft_control
3704 : TYPE(section_vals_type), POINTER :: hairy_probes_section
3705 :
3706 : INTEGER :: i, j, jj, kk, n_rep
3707 4 : INTEGER, DIMENSION(:), POINTER :: tmplist
3708 :
3709 12 : DO i = 1, SIZE(dft_control%probe)
3710 8 : NULLIFY (dft_control%probe(i)%atom_ids)
3711 :
3712 8 : CALL section_vals_val_get(hairy_probes_section, "ATOM_IDS", i_rep_section=i, n_rep_val=n_rep)
3713 8 : jj = 0
3714 16 : DO kk = 1, n_rep
3715 8 : CALL section_vals_val_get(hairy_probes_section, "ATOM_IDS", i_rep_section=i, i_rep_val=kk, i_vals=tmplist)
3716 16 : jj = jj + SIZE(tmplist)
3717 : END DO
3718 :
3719 8 : dft_control%probe(i)%natoms = jj
3720 8 : IF (dft_control%probe(i)%natoms < 1) THEN
3721 0 : CPABORT("Need at least 1 atom to use hair probes formalism")
3722 : END IF
3723 24 : ALLOCATE (dft_control%probe(i)%atom_ids(dft_control%probe(i)%natoms))
3724 :
3725 8 : jj = 0
3726 16 : DO kk = 1, n_rep
3727 8 : CALL section_vals_val_get(hairy_probes_section, "ATOM_IDS", i_rep_section=i, i_rep_val=kk, i_vals=tmplist)
3728 24 : DO j = 1, SIZE(tmplist)
3729 8 : jj = jj + 1
3730 16 : dft_control%probe(i)%atom_ids(jj) = tmplist(j)
3731 : END DO
3732 : END DO
3733 :
3734 8 : CALL section_vals_val_get(hairy_probes_section, "MU", i_rep_section=i, r_val=dft_control%probe(i)%mu)
3735 :
3736 8 : CALL section_vals_val_get(hairy_probes_section, "T", i_rep_section=i, r_val=dft_control%probe(i)%T)
3737 :
3738 8 : CALL section_vals_val_get(hairy_probes_section, "ALPHA", i_rep_section=i, r_val=dft_control%probe(i)%alpha)
3739 :
3740 20 : CALL section_vals_val_get(hairy_probes_section, "eps_hp", i_rep_section=i, r_val=dft_control%probe(i)%eps_hp)
3741 : END DO
3742 :
3743 4 : END SUBROUTINE read_hairy_probes_sections
3744 : ! **************************************************************************************************
3745 :
3746 : END MODULE cp_control_utils
|
