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