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