Line data Source code
1 : !--------------------------------------------------------------------------------------------------!
2 : ! CP2K: A general program to perform molecular dynamics simulations !
3 : ! Copyright 2000-2024 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, Dewar1977, Dewar1985, Elstner1998, Fattebert2002, Grimme2017, Hu2007, &
14 : Krack2000, Lippert1997, Lippert1999, Porezag1995, Repasky2002, Rocha2006, Schenter2008, &
15 : Seifert1996, Souza2002, Stengel2009, Stewart1989, Stewart2007, Thiel1992, Umari2002, &
16 : VanVoorhis2015, VandeVondele2005a, VandeVondele2005b, Yin2017, Zhechkov2005, cite_reference
17 : USE cp_control_types, ONLY: &
18 : admm_control_create, admm_control_type, ddapc_control_create, ddapc_restraint_type, &
19 : dft_control_create, dft_control_type, efield_type, expot_control_create, &
20 : maxwell_control_create, qs_control_type, tddfpt2_control_type, tddfpt_control_create, &
21 : tddfpt_control_type, xtb_control_type
22 : USE cp_files, ONLY: close_file,&
23 : open_file
24 : USE cp_log_handling, ONLY: cp_get_default_logger,&
25 : cp_logger_type
26 : USE cp_output_handling, ONLY: cp_print_key_finished_output,&
27 : cp_print_key_unit_nr
28 : USE cp_units, ONLY: cp_unit_from_cp2k,&
29 : cp_unit_to_cp2k
30 : USE force_fields_input, ONLY: read_gp_section
31 : USE input_constants, ONLY: &
32 : admm1_type, admm2_type, admmp_type, admmq_type, admms_type, constant_env, custom_env, &
33 : do_admm_aux_exch_func_bee, do_admm_aux_exch_func_bee_libxc, do_admm_aux_exch_func_default, &
34 : do_admm_aux_exch_func_default_libxc, do_admm_aux_exch_func_none, &
35 : do_admm_aux_exch_func_opt, do_admm_aux_exch_func_opt_libxc, do_admm_aux_exch_func_pbex, &
36 : do_admm_aux_exch_func_pbex_libxc, do_admm_aux_exch_func_sx_libxc, &
37 : do_admm_basis_projection, do_admm_blocked_projection, do_admm_blocking_purify_full, &
38 : do_admm_charge_constrained_projection, do_admm_exch_scaling_merlot, &
39 : do_admm_exch_scaling_none, do_admm_purify_cauchy, do_admm_purify_cauchy_subspace, &
40 : do_admm_purify_mcweeny, do_admm_purify_mo_diag, do_admm_purify_mo_no_diag, &
41 : do_admm_purify_none, do_admm_purify_none_dm, do_ddapc_constraint, do_ddapc_restraint, &
42 : do_method_am1, do_method_dftb, do_method_gapw, do_method_gapw_xc, do_method_gpw, &
43 : do_method_lrigpw, do_method_mndo, do_method_mndod, do_method_ofgpw, do_method_pdg, &
44 : do_method_pm3, do_method_pm6, do_method_pm6fm, do_method_pnnl, do_method_rigpw, &
45 : do_method_rm1, do_method_xtb, do_pwgrid_ns_fullspace, do_pwgrid_ns_halfspace, &
46 : do_pwgrid_spherical, do_s2_constraint, do_s2_restraint, do_se_is_kdso, do_se_is_kdso_d, &
47 : do_se_is_slater, do_se_lr_ewald, do_se_lr_ewald_gks, do_se_lr_ewald_r3, do_se_lr_none, &
48 : gapw_1c_large, gapw_1c_medium, gapw_1c_orb, gapw_1c_small, gapw_1c_very_large, &
49 : gaussian_env, no_admm_type, numerical, ramp_env, real_time_propagation, sccs_andreussi, &
50 : sccs_derivative_cd3, sccs_derivative_cd5, sccs_derivative_cd7, sccs_derivative_fft, &
51 : sccs_fattebert_gygi, sic_ad, sic_eo, sic_list_all, sic_list_unpaired, sic_mauri_spz, &
52 : sic_mauri_us, sic_none, slater, tddfpt_dipole_length, tddfpt_excitations, &
53 : tddfpt_kernel_stda, use_mom_ref_user
54 : USE input_cp2k_check, ONLY: xc_functionals_expand
55 : USE input_cp2k_dft, ONLY: create_dft_section
56 : USE input_enumeration_types, ONLY: enum_i2c,&
57 : enumeration_type
58 : USE input_keyword_types, ONLY: keyword_get,&
59 : keyword_type
60 : USE input_section_types, ONLY: &
61 : section_get_ival, section_get_keyword, section_release, section_type, section_vals_get, &
62 : section_vals_get_subs_vals, section_vals_type, section_vals_val_get, section_vals_val_set
63 : USE kinds, ONLY: default_path_length,&
64 : default_string_length,&
65 : dp
66 : USE mathconstants, ONLY: fourpi
67 : USE pair_potential_types, ONLY: pair_potential_reallocate
68 : USE periodic_table, ONLY: get_ptable_info
69 : USE qs_cdft_utils, ONLY: read_cdft_control_section
70 : USE string_utilities, ONLY: uppercase
71 : USE util, ONLY: sort
72 : USE xc, ONLY: xc_uses_kinetic_energy_density,&
73 : xc_uses_norm_drho
74 : USE xc_input_constants, ONLY: xc_deriv_collocate
75 : USE xc_write_output, ONLY: xc_write
76 : #include "./base/base_uses.f90"
77 :
78 : IMPLICIT NONE
79 :
80 : PRIVATE
81 :
82 : CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'cp_control_utils'
83 :
84 : PUBLIC :: read_dft_control, &
85 : read_mgrid_section, &
86 : read_qs_section, &
87 : read_tddfpt_control, &
88 : read_tddfpt2_control, &
89 : write_dft_control, &
90 : write_qs_control, &
91 : write_admm_control, &
92 : read_ddapc_section
93 : CONTAINS
94 :
95 : ! **************************************************************************************************
96 : !> \brief ...
97 : !> \param dft_control ...
98 : !> \param dft_section ...
99 : ! **************************************************************************************************
100 78600 : SUBROUTINE read_dft_control(dft_control, dft_section)
101 : TYPE(dft_control_type), POINTER :: dft_control
102 : TYPE(section_vals_type), POINTER :: dft_section
103 :
104 : CHARACTER(len=default_path_length) :: basis_set_file_name, potential_file_name
105 : CHARACTER(LEN=default_string_length), &
106 6550 : DIMENSION(:), POINTER :: tmpstringlist
107 : INTEGER :: admmtype, excitations, irep, isize, &
108 : method_id, nrep, xc_deriv_method_id
109 : LOGICAL :: do_ot, do_rtp, exopt1, exopt2, exopt3, &
110 : explicit, is_present, l_param, not_SE, &
111 : was_present
112 : REAL(KIND=dp) :: density_cut, gradient_cut, tau_cut
113 6550 : REAL(KIND=dp), DIMENSION(:), POINTER :: pol
114 : TYPE(cp_logger_type), POINTER :: logger
115 : TYPE(section_vals_type), POINTER :: maxwell_section, sccs_section, &
116 : scf_section, tmp_section, &
117 : xc_fun_section, xc_section
118 :
119 6550 : was_present = .FALSE.
120 :
121 6550 : logger => cp_get_default_logger()
122 :
123 6550 : NULLIFY (tmp_section, xc_fun_section, xc_section)
124 6550 : ALLOCATE (dft_control)
125 6550 : CALL dft_control_create(dft_control)
126 : ! determine wheather this is a semiempirical or DFTB run
127 : ! --> (no XC section needs to be provided)
128 6550 : not_SE = .TRUE.
129 6550 : CALL section_vals_val_get(dft_section, "QS%METHOD", i_val=method_id)
130 1468 : SELECT CASE (method_id)
131 : CASE (do_method_dftb, do_method_xtb, do_method_mndo, do_method_am1, do_method_pm3, do_method_pnnl, &
132 : do_method_pm6, do_method_pm6fm, do_method_pdg, do_method_rm1, do_method_mndod)
133 6550 : not_SE = .FALSE.
134 : END SELECT
135 : ! Check for XC section and XC_FUNCTIONAL section
136 6550 : xc_section => section_vals_get_subs_vals(dft_section, "XC")
137 6550 : CALL section_vals_get(xc_section, explicit=is_present)
138 6550 : IF (.NOT. is_present .AND. not_SE) THEN
139 0 : CPABORT("XC section missing.")
140 : END IF
141 6550 : IF (is_present) THEN
142 5094 : CALL section_vals_val_get(xc_section, "density_cutoff", r_val=density_cut)
143 5094 : CALL section_vals_val_get(xc_section, "gradient_cutoff", r_val=gradient_cut)
144 5094 : CALL section_vals_val_get(xc_section, "tau_cutoff", r_val=tau_cut)
145 : ! Perform numerical stability checks and possibly correct the issues
146 5094 : IF (density_cut <= EPSILON(0.0_dp)*100.0_dp) &
147 : CALL cp_warn(__LOCATION__, &
148 : "DENSITY_CUTOFF lower than 100*EPSILON, where EPSILON is the machine precision. "// &
149 0 : "This may lead to numerical problems. Setting up shake_tol to 100*EPSILON! ")
150 5094 : density_cut = MAX(EPSILON(0.0_dp)*100.0_dp, density_cut)
151 5094 : IF (gradient_cut <= EPSILON(0.0_dp)*100.0_dp) &
152 : CALL cp_warn(__LOCATION__, &
153 : "GRADIENT_CUTOFF lower than 100*EPSILON, where EPSILON is the machine precision. "// &
154 0 : "This may lead to numerical problems. Setting up shake_tol to 100*EPSILON! ")
155 5094 : gradient_cut = MAX(EPSILON(0.0_dp)*100.0_dp, gradient_cut)
156 5094 : IF (tau_cut <= EPSILON(0.0_dp)*100.0_dp) &
157 : CALL cp_warn(__LOCATION__, &
158 : "TAU_CUTOFF lower than 100*EPSILON, where EPSILON is the machine precision. "// &
159 0 : "This may lead to numerical problems. Setting up shake_tol to 100*EPSILON! ")
160 5094 : tau_cut = MAX(EPSILON(0.0_dp)*100.0_dp, tau_cut)
161 5094 : CALL section_vals_val_set(xc_section, "density_cutoff", r_val=density_cut)
162 5094 : CALL section_vals_val_set(xc_section, "gradient_cutoff", r_val=gradient_cut)
163 5094 : CALL section_vals_val_set(xc_section, "tau_cutoff", r_val=tau_cut)
164 : END IF
165 6550 : xc_fun_section => section_vals_get_subs_vals(xc_section, "XC_FUNCTIONAL")
166 6550 : CALL section_vals_get(xc_fun_section, explicit=is_present)
167 6550 : IF (.NOT. is_present .AND. not_SE) THEN
168 0 : CPABORT("XC_FUNCTIONAL section missing.")
169 : END IF
170 6550 : scf_section => section_vals_get_subs_vals(dft_section, "SCF")
171 6550 : CALL section_vals_val_get(dft_section, "UKS", l_val=dft_control%uks)
172 6550 : CALL section_vals_val_get(dft_section, "ROKS", l_val=dft_control%roks)
173 6550 : IF (dft_control%uks .OR. dft_control%roks) THEN
174 1471 : dft_control%nspins = 2
175 : ELSE
176 5079 : dft_control%nspins = 1
177 : END IF
178 :
179 6550 : dft_control%lsd = (dft_control%nspins > 1)
180 6550 : dft_control%use_kinetic_energy_density = xc_uses_kinetic_energy_density(xc_fun_section, dft_control%lsd)
181 :
182 6550 : xc_deriv_method_id = section_get_ival(xc_section, "XC_GRID%XC_DERIV")
183 : dft_control%drho_by_collocation = (xc_uses_norm_drho(xc_fun_section, dft_control%lsd) &
184 6550 : .AND. (xc_deriv_method_id == xc_deriv_collocate))
185 6550 : IF (dft_control%drho_by_collocation) THEN
186 0 : CPABORT("derivatives by collocation not implemented")
187 : END IF
188 :
189 : ! Automatic auxiliary basis set generation
190 6550 : CALL section_vals_val_get(dft_section, "AUTO_BASIS", n_rep_val=nrep)
191 13100 : DO irep = 1, nrep
192 6550 : CALL section_vals_val_get(dft_section, "AUTO_BASIS", i_rep_val=irep, c_vals=tmpstringlist)
193 13100 : IF (SIZE(tmpstringlist) == 2) THEN
194 6550 : CALL uppercase(tmpstringlist(2))
195 6550 : SELECT CASE (tmpstringlist(2))
196 : CASE ("X")
197 78 : isize = -1
198 : CASE ("SMALL")
199 78 : isize = 0
200 : CASE ("MEDIUM")
201 48 : isize = 1
202 : CASE ("LARGE")
203 0 : isize = 2
204 : CASE ("HUGE")
205 2 : isize = 3
206 : CASE DEFAULT
207 6550 : CPWARN("Unknown basis size in AUTO_BASIS keyword:"//TRIM(tmpstringlist(1)))
208 : END SELECT
209 : !
210 6552 : SELECT CASE (tmpstringlist(1))
211 : CASE ("X")
212 : CASE ("RI_AUX")
213 2 : dft_control%auto_basis_ri_aux = isize
214 : CASE ("AUX_FIT")
215 0 : dft_control%auto_basis_aux_fit = isize
216 : CASE ("LRI_AUX")
217 0 : dft_control%auto_basis_lri_aux = isize
218 : CASE ("P_LRI_AUX")
219 0 : dft_control%auto_basis_p_lri_aux = isize
220 : CASE ("RI_HXC")
221 0 : dft_control%auto_basis_ri_hxc = isize
222 : CASE ("RI_XAS")
223 48 : dft_control%auto_basis_ri_xas = isize
224 : CASE ("RI_HFX")
225 78 : dft_control%auto_basis_ri_hfx = isize
226 : CASE DEFAULT
227 6550 : CPWARN("Unknown basis type in AUTO_BASIS keyword:"//TRIM(tmpstringlist(1)))
228 : END SELECT
229 : ELSE
230 : CALL cp_abort(__LOCATION__, &
231 0 : "AUTO_BASIS keyword in &DFT section has a wrong number of arguments.")
232 : END IF
233 : END DO
234 :
235 : !! check if we do wavefunction fitting
236 6550 : tmp_section => section_vals_get_subs_vals(dft_section, "AUXILIARY_DENSITY_MATRIX_METHOD")
237 6550 : CALL section_vals_get(tmp_section, explicit=is_present)
238 6550 : dft_control%do_admm = is_present
239 6550 : dft_control%do_admm_mo = .FALSE.
240 6550 : dft_control%do_admm_dm = .FALSE.
241 6550 : IF (is_present) THEN
242 : do_ot = .FALSE.
243 442 : CALL section_vals_val_get(scf_section, "OT%_SECTION_PARAMETERS_", l_val=do_ot)
244 442 : CALL admm_control_create(dft_control%admm_control)
245 :
246 442 : CALL section_vals_val_get(dft_section, "AUXILIARY_DENSITY_MATRIX_METHOD%ADMM_TYPE", i_val=admmtype)
247 442 : CALL section_vals_val_get(dft_section, "AUXILIARY_DENSITY_MATRIX_METHOD%ADMM_PURIFICATION_METHOD", explicit=exopt1)
248 442 : CALL section_vals_val_get(dft_section, "AUXILIARY_DENSITY_MATRIX_METHOD%METHOD", explicit=exopt2)
249 442 : CALL section_vals_val_get(dft_section, "AUXILIARY_DENSITY_MATRIX_METHOD%EXCH_SCALING_MODEL", explicit=exopt3)
250 442 : dft_control%admm_control%admm_type = admmtype
251 432 : SELECT CASE (admmtype)
252 : CASE (no_admm_type)
253 432 : CALL section_vals_val_get(dft_section, "AUXILIARY_DENSITY_MATRIX_METHOD%ADMM_PURIFICATION_METHOD", i_val=method_id)
254 432 : dft_control%admm_control%purification_method = method_id
255 432 : CALL section_vals_val_get(dft_section, "AUXILIARY_DENSITY_MATRIX_METHOD%METHOD", i_val=method_id)
256 432 : dft_control%admm_control%method = method_id
257 432 : CALL section_vals_val_get(dft_section, "AUXILIARY_DENSITY_MATRIX_METHOD%EXCH_SCALING_MODEL", i_val=method_id)
258 432 : dft_control%admm_control%scaling_model = method_id
259 : CASE (admm1_type)
260 : ! METHOD BASIS_PROJECTION
261 : ! ADMM_PURIFICATION_METHOD choose
262 : ! EXCH_SCALING_MODEL NONE
263 2 : CALL section_vals_val_get(dft_section, "AUXILIARY_DENSITY_MATRIX_METHOD%ADMM_PURIFICATION_METHOD", i_val=method_id)
264 2 : dft_control%admm_control%purification_method = method_id
265 2 : dft_control%admm_control%method = do_admm_basis_projection
266 2 : dft_control%admm_control%scaling_model = do_admm_exch_scaling_none
267 : CASE (admm2_type)
268 : ! METHOD BASIS_PROJECTION
269 : ! ADMM_PURIFICATION_METHOD NONE
270 : ! EXCH_SCALING_MODEL NONE
271 2 : dft_control%admm_control%purification_method = do_admm_purify_none
272 2 : dft_control%admm_control%method = do_admm_basis_projection
273 2 : dft_control%admm_control%scaling_model = do_admm_exch_scaling_none
274 : CASE (admms_type)
275 : ! ADMM_PURIFICATION_METHOD NONE
276 : ! METHOD CHARGE_CONSTRAINED_PROJECTION
277 : ! EXCH_SCALING_MODEL MERLOT
278 2 : dft_control%admm_control%purification_method = do_admm_purify_none
279 2 : dft_control%admm_control%method = do_admm_charge_constrained_projection
280 2 : dft_control%admm_control%scaling_model = do_admm_exch_scaling_merlot
281 : CASE (admmp_type)
282 : ! ADMM_PURIFICATION_METHOD NONE
283 : ! METHOD BASIS_PROJECTION
284 : ! EXCH_SCALING_MODEL MERLOT
285 2 : dft_control%admm_control%purification_method = do_admm_purify_none
286 2 : dft_control%admm_control%method = do_admm_basis_projection
287 2 : dft_control%admm_control%scaling_model = do_admm_exch_scaling_merlot
288 : CASE (admmq_type)
289 : ! ADMM_PURIFICATION_METHOD NONE
290 : ! METHOD CHARGE_CONSTRAINED_PROJECTION
291 : ! EXCH_SCALING_MODEL NONE
292 2 : dft_control%admm_control%purification_method = do_admm_purify_none
293 2 : dft_control%admm_control%method = do_admm_charge_constrained_projection
294 2 : dft_control%admm_control%scaling_model = do_admm_exch_scaling_none
295 : CASE DEFAULT
296 : CALL cp_abort(__LOCATION__, &
297 442 : "ADMM_TYPE keyword in &AUXILIARY_DENSITY_MATRIX_METHOD section has a wrong value.")
298 : END SELECT
299 :
300 : CALL section_vals_val_get(dft_section, "AUXILIARY_DENSITY_MATRIX_METHOD%EPS_FILTER", &
301 442 : r_val=dft_control%admm_control%eps_filter)
302 :
303 442 : CALL section_vals_val_get(dft_section, "AUXILIARY_DENSITY_MATRIX_METHOD%EXCH_CORRECTION_FUNC", i_val=method_id)
304 442 : dft_control%admm_control%aux_exch_func = method_id
305 :
306 : ! parameters for X functional
307 442 : dft_control%admm_control%aux_exch_func_param = .FALSE.
308 : CALL section_vals_val_get(dft_section, "AUXILIARY_DENSITY_MATRIX_METHOD%OPTX_A1", explicit=explicit, &
309 442 : r_val=dft_control%admm_control%aux_x_param(1))
310 442 : IF (explicit) dft_control%admm_control%aux_exch_func_param = .TRUE.
311 : CALL section_vals_val_get(dft_section, "AUXILIARY_DENSITY_MATRIX_METHOD%OPTX_A2", explicit=explicit, &
312 442 : r_val=dft_control%admm_control%aux_x_param(2))
313 442 : IF (explicit) dft_control%admm_control%aux_exch_func_param = .TRUE.
314 : CALL section_vals_val_get(dft_section, "AUXILIARY_DENSITY_MATRIX_METHOD%OPTX_GAMMA", explicit=explicit, &
315 442 : r_val=dft_control%admm_control%aux_x_param(3))
316 442 : IF (explicit) dft_control%admm_control%aux_exch_func_param = .TRUE.
317 :
318 442 : CALL read_admm_block_list(dft_control%admm_control, dft_section)
319 :
320 : ! check for double assignments
321 2 : SELECT CASE (admmtype)
322 : CASE (admm2_type)
323 2 : IF (exopt2) CALL cp_warn(__LOCATION__, &
324 0 : "Value of ADMM_PURIFICATION_METHOD keyword will be overwritten with ADMM_TYPE selections.")
325 2 : IF (exopt3) CALL cp_warn(__LOCATION__, &
326 0 : "Value of EXCH_SCALING_MODEL keyword will be overwritten with ADMM_TYPE selections.")
327 : CASE (admm1_type, admms_type, admmp_type, admmq_type)
328 8 : IF (exopt1) CALL cp_warn(__LOCATION__, &
329 0 : "Value of METHOD keyword will be overwritten with ADMM_TYPE selections.")
330 8 : IF (exopt2) CALL cp_warn(__LOCATION__, &
331 0 : "Value of METHOD keyword will be overwritten with ADMM_TYPE selections.")
332 8 : IF (exopt3) CALL cp_warn(__LOCATION__, &
333 442 : "Value of EXCH_SCALING_MODEL keyword will be overwritten with ADMM_TYPE selections.")
334 : END SELECT
335 :
336 : ! In the case of charge-constrained projection (e.g. according to Merlot),
337 : ! there is no purification needed and hence, do_admm_purify_none has to be set.
338 :
339 : IF ((dft_control%admm_control%method == do_admm_blocking_purify_full .OR. &
340 : dft_control%admm_control%method == do_admm_blocked_projection) &
341 442 : .AND. dft_control%admm_control%scaling_model == do_admm_exch_scaling_merlot) THEN
342 0 : CPABORT("ADMM: Blocking and Merlot scaling are mutually exclusive.")
343 : END IF
344 :
345 442 : IF (dft_control%admm_control%method == do_admm_charge_constrained_projection .AND. &
346 : dft_control%admm_control%purification_method /= do_admm_purify_none) THEN
347 : CALL cp_abort(__LOCATION__, &
348 : "ADMM: In the case of METHOD=CHARGE_CONSTRAINED_PROJECTION, "// &
349 0 : "ADMM_PURIFICATION_METHOD=NONE has to be set.")
350 : END IF
351 :
352 442 : IF (dft_control%admm_control%purification_method == do_admm_purify_mo_diag .OR. &
353 : dft_control%admm_control%purification_method == do_admm_purify_mo_no_diag) THEN
354 60 : IF (dft_control%admm_control%method /= do_admm_basis_projection) &
355 0 : CPABORT("ADMM: Chosen purification requires BASIS_PROJECTION")
356 :
357 60 : IF (.NOT. do_ot) CPABORT("ADMM: MO-based purification requires OT.")
358 : END IF
359 :
360 442 : IF (dft_control%admm_control%purification_method == do_admm_purify_none_dm .OR. &
361 : dft_control%admm_control%purification_method == do_admm_purify_mcweeny) THEN
362 14 : dft_control%do_admm_dm = .TRUE.
363 : ELSE
364 428 : dft_control%do_admm_mo = .TRUE.
365 : END IF
366 : END IF
367 :
368 : ! Set restricted to true, if both OT and ROKS are requested
369 : !MK in principle dft_control%restricted could be dropped completely like the
370 : !MK input key by using only dft_control%roks now
371 6550 : CALL section_vals_val_get(scf_section, "OT%_SECTION_PARAMETERS_", l_val=l_param)
372 6550 : dft_control%restricted = (dft_control%roks .AND. l_param)
373 :
374 6550 : CALL section_vals_val_get(dft_section, "CHARGE", i_val=dft_control%charge)
375 6550 : CALL section_vals_val_get(dft_section, "MULTIPLICITY", i_val=dft_control%multiplicity)
376 6550 : CALL section_vals_val_get(dft_section, "RELAX_MULTIPLICITY", r_val=dft_control%relax_multiplicity)
377 6550 : IF (dft_control%relax_multiplicity > 0.0_dp) THEN
378 8 : IF (.NOT. dft_control%uks) &
379 : CALL cp_abort(__LOCATION__, "The option RELAX_MULTIPLICITY is only valid for "// &
380 0 : "unrestricted Kohn-Sham (UKS) calculations")
381 : END IF
382 :
383 : ! check for the presence of the low spin roks section
384 6550 : tmp_section => section_vals_get_subs_vals(dft_section, "LOW_SPIN_ROKS")
385 6550 : CALL section_vals_get(tmp_section, explicit=dft_control%low_spin_roks)
386 :
387 6550 : dft_control%sic_method_id = sic_none
388 6550 : dft_control%sic_scaling_a = 1.0_dp
389 6550 : dft_control%sic_scaling_b = 1.0_dp
390 :
391 : ! DFT+U
392 6550 : dft_control%dft_plus_u = .FALSE.
393 6550 : CALL section_vals_val_get(dft_section, "PLUS_U_METHOD", i_val=method_id)
394 6550 : dft_control%plus_u_method_id = method_id
395 :
396 : ! Smearing in use
397 6550 : dft_control%smear = .FALSE.
398 :
399 : ! Surface dipole correction
400 6550 : dft_control%correct_surf_dip = .FALSE.
401 6550 : CALL section_vals_val_get(dft_section, "SURFACE_DIPOLE_CORRECTION", l_val=dft_control%correct_surf_dip)
402 6550 : CALL section_vals_val_get(dft_section, "SURF_DIP_DIR", i_val=dft_control%dir_surf_dip)
403 6550 : dft_control%pos_dir_surf_dip = -1.0_dp
404 6550 : CALL section_vals_val_get(dft_section, "SURF_DIP_POS", r_val=dft_control%pos_dir_surf_dip)
405 : ! another logical variable, surf_dip_correct_switch, is introduced for
406 : ! implementation of "SURF_DIP_SWITCH" [SGh]
407 6550 : dft_control%switch_surf_dip = .FALSE.
408 6550 : dft_control%surf_dip_correct_switch = dft_control%correct_surf_dip
409 6550 : CALL section_vals_val_get(dft_section, "SURF_DIP_SWITCH", l_val=dft_control%switch_surf_dip)
410 6550 : dft_control%correct_el_density_dip = .FALSE.
411 6550 : CALL section_vals_val_get(dft_section, "CORE_CORR_DIP", l_val=dft_control%correct_el_density_dip)
412 6550 : IF (dft_control%correct_el_density_dip) THEN
413 4 : IF (dft_control%correct_surf_dip) THEN
414 : ! Do nothing, move on
415 : ELSE
416 0 : dft_control%correct_el_density_dip = .FALSE.
417 0 : CPWARN("CORE_CORR_DIP keyword is activated only if SURFACE_DIPOLE_CORRECTION is TRUE")
418 : END IF
419 : END IF
420 :
421 : CALL section_vals_val_get(dft_section, "BASIS_SET_FILE_NAME", &
422 6550 : c_val=basis_set_file_name)
423 : CALL section_vals_val_get(dft_section, "POTENTIAL_FILE_NAME", &
424 6550 : c_val=potential_file_name)
425 :
426 : ! Read the input section
427 6550 : tmp_section => section_vals_get_subs_vals(dft_section, "sic")
428 : CALL section_vals_val_get(tmp_section, "SIC_METHOD", &
429 6550 : i_val=dft_control%sic_method_id)
430 : CALL section_vals_val_get(tmp_section, "ORBITAL_SET", &
431 6550 : i_val=dft_control%sic_list_id)
432 : CALL section_vals_val_get(tmp_section, "SIC_SCALING_A", &
433 6550 : r_val=dft_control%sic_scaling_a)
434 : CALL section_vals_val_get(tmp_section, "SIC_SCALING_B", &
435 6550 : r_val=dft_control%sic_scaling_b)
436 :
437 : ! Determine if this is a TDDFPT run
438 6550 : CALL section_vals_val_get(dft_section, "EXCITATIONS", i_val=excitations)
439 6550 : dft_control%do_tddfpt_calculation = (excitations == tddfpt_excitations)
440 6550 : IF (dft_control%do_tddfpt_calculation) THEN
441 12 : CALL tddfpt_control_create(dft_control%tddfpt_control)
442 : END IF
443 :
444 6550 : do_rtp = .FALSE.
445 6550 : tmp_section => section_vals_get_subs_vals(dft_section, "REAL_TIME_PROPAGATION")
446 6550 : CALL section_vals_get(tmp_section, explicit=is_present)
447 6550 : IF (is_present) THEN
448 236 : CALL read_rtp_section(dft_control, tmp_section)
449 236 : do_rtp = .TRUE.
450 : END IF
451 :
452 : ! Read the input section
453 6550 : tmp_section => section_vals_get_subs_vals(dft_section, "XAS")
454 6550 : CALL section_vals_get(tmp_section, explicit=dft_control%do_xas_calculation)
455 6550 : IF (dft_control%do_xas_calculation) THEN
456 : ! Override with section parameter
457 : CALL section_vals_val_get(tmp_section, "_SECTION_PARAMETERS_", &
458 42 : l_val=dft_control%do_xas_calculation)
459 : END IF
460 :
461 6550 : tmp_section => section_vals_get_subs_vals(dft_section, "XAS_TDP")
462 6550 : CALL section_vals_get(tmp_section, explicit=dft_control%do_xas_tdp_calculation)
463 6550 : IF (dft_control%do_xas_tdp_calculation) THEN
464 : ! Override with section parameter
465 : CALL section_vals_val_get(tmp_section, "_SECTION_PARAMETERS_", &
466 50 : l_val=dft_control%do_xas_tdp_calculation)
467 : END IF
468 :
469 : ! Read the finite field input section
470 6550 : dft_control%apply_efield = .FALSE.
471 6550 : dft_control%apply_efield_field = .FALSE. !this is for RTP
472 6550 : dft_control%apply_vector_potential = .FALSE. !this is for RTP
473 6550 : tmp_section => section_vals_get_subs_vals(dft_section, "EFIELD")
474 6550 : CALL section_vals_get(tmp_section, n_repetition=nrep, explicit=is_present)
475 6550 : IF (is_present) THEN
476 976 : ALLOCATE (dft_control%efield_fields(nrep))
477 244 : CALL read_efield_sections(dft_control, tmp_section)
478 244 : IF (do_rtp) THEN
479 18 : IF (.NOT. dft_control%rtp_control%velocity_gauge) THEN
480 10 : dft_control%apply_efield_field = .TRUE.
481 : ELSE
482 8 : dft_control%apply_vector_potential = .TRUE.
483 : ! Use this input value of vector potential to (re)start RTP
484 32 : dft_control%rtp_control%vec_pot = dft_control%efield_fields(1)%efield%vec_pot_initial
485 : END IF
486 : ELSE
487 226 : dft_control%apply_efield = .TRUE.
488 226 : CPASSERT(nrep == 1)
489 : END IF
490 : END IF
491 :
492 : ! Read the finite field input section for periodic fields
493 6550 : tmp_section => section_vals_get_subs_vals(dft_section, "PERIODIC_EFIELD")
494 6550 : CALL section_vals_get(tmp_section, explicit=dft_control%apply_period_efield)
495 6550 : IF (dft_control%apply_period_efield) THEN
496 392 : ALLOCATE (dft_control%period_efield)
497 56 : CALL section_vals_val_get(tmp_section, "POLARISATION", r_vals=pol)
498 392 : dft_control%period_efield%polarisation(1:3) = pol(1:3)
499 56 : CALL section_vals_val_get(tmp_section, "D_FILTER", r_vals=pol)
500 392 : dft_control%period_efield%d_filter(1:3) = pol(1:3)
501 : CALL section_vals_val_get(tmp_section, "INTENSITY", &
502 56 : r_val=dft_control%period_efield%strength)
503 56 : dft_control%period_efield%displacement_field = .FALSE.
504 : CALL section_vals_val_get(tmp_section, "DISPLACEMENT_FIELD", &
505 56 : l_val=dft_control%period_efield%displacement_field)
506 : ! periodic fields don't work with RTP
507 56 : CPASSERT(.NOT. do_rtp)
508 56 : IF (dft_control%period_efield%displacement_field) THEN
509 16 : CALL cite_reference(Stengel2009)
510 : ELSE
511 40 : CALL cite_reference(Souza2002)
512 40 : CALL cite_reference(Umari2002)
513 : END IF
514 : END IF
515 :
516 : ! Read the external potential input section
517 6550 : tmp_section => section_vals_get_subs_vals(dft_section, "EXTERNAL_POTENTIAL")
518 6550 : CALL section_vals_get(tmp_section, explicit=dft_control%apply_external_potential)
519 6550 : IF (dft_control%apply_external_potential) THEN
520 16 : CALL expot_control_create(dft_control%expot_control)
521 : CALL section_vals_val_get(tmp_section, "READ_FROM_CUBE", &
522 16 : l_val=dft_control%expot_control%read_from_cube)
523 : CALL section_vals_val_get(tmp_section, "STATIC", &
524 16 : l_val=dft_control%expot_control%static)
525 : CALL section_vals_val_get(tmp_section, "SCALING_FACTOR", &
526 16 : r_val=dft_control%expot_control%scaling_factor)
527 : ! External potential using Maxwell equation
528 16 : maxwell_section => section_vals_get_subs_vals(tmp_section, "MAXWELL")
529 16 : CALL section_vals_get(maxwell_section, explicit=is_present)
530 16 : IF (is_present) THEN
531 0 : dft_control%expot_control%maxwell_solver = .TRUE.
532 0 : CALL maxwell_control_create(dft_control%maxwell_control)
533 : ! read the input values from Maxwell section
534 : CALL section_vals_val_get(maxwell_section, "TEST_REAL", &
535 0 : r_val=dft_control%maxwell_control%real_test)
536 : CALL section_vals_val_get(maxwell_section, "TEST_INTEGER", &
537 0 : i_val=dft_control%maxwell_control%int_test)
538 : CALL section_vals_val_get(maxwell_section, "TEST_LOGICAL", &
539 0 : l_val=dft_control%maxwell_control%log_test)
540 : ELSE
541 16 : dft_control%expot_control%maxwell_solver = .FALSE.
542 : END IF
543 : END IF
544 :
545 : ! Read the SCCS input section if present
546 6550 : sccs_section => section_vals_get_subs_vals(dft_section, "SCCS")
547 6550 : CALL section_vals_get(sccs_section, explicit=is_present)
548 6550 : IF (is_present) THEN
549 : ! Check section parameter if SCCS is activated
550 : CALL section_vals_val_get(sccs_section, "_SECTION_PARAMETERS_", &
551 10 : l_val=dft_control%do_sccs)
552 10 : IF (dft_control%do_sccs) THEN
553 10 : ALLOCATE (dft_control%sccs_control)
554 : CALL section_vals_val_get(sccs_section, "RELATIVE_PERMITTIVITY", &
555 10 : r_val=dft_control%sccs_control%epsilon_solvent)
556 : CALL section_vals_val_get(sccs_section, "ALPHA", &
557 10 : r_val=dft_control%sccs_control%alpha_solvent)
558 : CALL section_vals_val_get(sccs_section, "BETA", &
559 10 : r_val=dft_control%sccs_control%beta_solvent)
560 : CALL section_vals_val_get(sccs_section, "DELTA_RHO", &
561 10 : r_val=dft_control%sccs_control%delta_rho)
562 : CALL section_vals_val_get(sccs_section, "DERIVATIVE_METHOD", &
563 10 : i_val=dft_control%sccs_control%derivative_method)
564 : CALL section_vals_val_get(sccs_section, "METHOD", &
565 10 : i_val=dft_control%sccs_control%method_id)
566 : CALL section_vals_val_get(sccs_section, "EPS_SCCS", &
567 10 : r_val=dft_control%sccs_control%eps_sccs)
568 : CALL section_vals_val_get(sccs_section, "EPS_SCF", &
569 10 : r_val=dft_control%sccs_control%eps_scf)
570 : CALL section_vals_val_get(sccs_section, "GAMMA", &
571 10 : r_val=dft_control%sccs_control%gamma_solvent)
572 : CALL section_vals_val_get(sccs_section, "MAX_ITER", &
573 10 : i_val=dft_control%sccs_control%max_iter)
574 : CALL section_vals_val_get(sccs_section, "MIXING", &
575 10 : r_val=dft_control%sccs_control%mixing)
576 18 : SELECT CASE (dft_control%sccs_control%method_id)
577 : CASE (sccs_andreussi)
578 8 : tmp_section => section_vals_get_subs_vals(sccs_section, "ANDREUSSI")
579 : CALL section_vals_val_get(tmp_section, "RHO_MAX", &
580 8 : r_val=dft_control%sccs_control%rho_max)
581 : CALL section_vals_val_get(tmp_section, "RHO_MIN", &
582 8 : r_val=dft_control%sccs_control%rho_min)
583 8 : IF (dft_control%sccs_control%rho_max < dft_control%sccs_control%rho_min) THEN
584 : CALL cp_abort(__LOCATION__, &
585 : "The SCCS parameter RHO_MAX is smaller than RHO_MIN. "// &
586 0 : "Please, check your input!")
587 : END IF
588 8 : CALL cite_reference(Andreussi2012)
589 : CASE (sccs_fattebert_gygi)
590 2 : tmp_section => section_vals_get_subs_vals(sccs_section, "FATTEBERT-GYGI")
591 : CALL section_vals_val_get(tmp_section, "BETA", &
592 2 : r_val=dft_control%sccs_control%beta)
593 2 : IF (dft_control%sccs_control%beta < 0.5_dp) THEN
594 : CALL cp_abort(__LOCATION__, &
595 : "A value smaller than 0.5 for the SCCS parameter beta "// &
596 0 : "causes numerical problems. Please, check your input!")
597 : END IF
598 : CALL section_vals_val_get(tmp_section, "RHO_ZERO", &
599 2 : r_val=dft_control%sccs_control%rho_zero)
600 2 : CALL cite_reference(Fattebert2002)
601 : CASE DEFAULT
602 10 : CPABORT("Invalid SCCS model specified. Please, check your input!")
603 : END SELECT
604 10 : CALL cite_reference(Yin2017)
605 : END IF
606 : END IF
607 :
608 : ! ZMP added input sections
609 : ! Read the external density input section
610 6550 : tmp_section => section_vals_get_subs_vals(dft_section, "EXTERNAL_DENSITY")
611 6550 : CALL section_vals_get(tmp_section, explicit=dft_control%apply_external_density)
612 :
613 : ! Read the external vxc input section
614 6550 : tmp_section => section_vals_get_subs_vals(dft_section, "EXTERNAL_VXC")
615 6550 : CALL section_vals_get(tmp_section, explicit=dft_control%apply_external_vxc)
616 :
617 6550 : END SUBROUTINE read_dft_control
618 :
619 : ! **************************************************************************************************
620 : !> \brief ...
621 : !> \param qs_control ...
622 : !> \param dft_section ...
623 : ! **************************************************************************************************
624 6550 : SUBROUTINE read_mgrid_section(qs_control, dft_section)
625 :
626 : TYPE(qs_control_type), INTENT(INOUT) :: qs_control
627 : TYPE(section_vals_type), POINTER :: dft_section
628 :
629 : CHARACTER(len=*), PARAMETER :: routineN = 'read_mgrid_section'
630 :
631 : INTEGER :: handle, igrid_level, ngrid_level
632 : LOGICAL :: explicit, multigrid_set
633 : REAL(dp) :: cutoff
634 6550 : REAL(dp), DIMENSION(:), POINTER :: cutofflist
635 : TYPE(section_vals_type), POINTER :: mgrid_section
636 :
637 6550 : CALL timeset(routineN, handle)
638 :
639 6550 : NULLIFY (mgrid_section, cutofflist)
640 6550 : mgrid_section => section_vals_get_subs_vals(dft_section, "MGRID")
641 :
642 6550 : CALL section_vals_val_get(mgrid_section, "NGRIDS", i_val=ngrid_level)
643 6550 : CALL section_vals_val_get(mgrid_section, "MULTIGRID_SET", l_val=multigrid_set)
644 6550 : CALL section_vals_val_get(mgrid_section, "CUTOFF", r_val=cutoff)
645 6550 : CALL section_vals_val_get(mgrid_section, "PROGRESSION_FACTOR", r_val=qs_control%progression_factor)
646 6550 : CALL section_vals_val_get(mgrid_section, "COMMENSURATE", l_val=qs_control%commensurate_mgrids)
647 6550 : CALL section_vals_val_get(mgrid_section, "REALSPACE", l_val=qs_control%realspace_mgrids)
648 6550 : CALL section_vals_val_get(mgrid_section, "REL_CUTOFF", r_val=qs_control%relative_cutoff)
649 : CALL section_vals_val_get(mgrid_section, "SKIP_LOAD_BALANCE_DISTRIBUTED", &
650 6550 : l_val=qs_control%skip_load_balance_distributed)
651 :
652 : ! For SE and DFTB possibly override with new defaults
653 6550 : IF (qs_control%semi_empirical .OR. qs_control%dftb .OR. qs_control%xtb) THEN
654 1468 : ngrid_level = 1
655 1468 : multigrid_set = .FALSE.
656 : ! Override default cutoff value unless user specified an explicit argument..
657 1468 : CALL section_vals_val_get(mgrid_section, "CUTOFF", explicit=explicit, r_val=cutoff)
658 1468 : IF (.NOT. explicit) cutoff = 1.0_dp
659 : END IF
660 :
661 19650 : ALLOCATE (qs_control%e_cutoff(ngrid_level))
662 6550 : qs_control%cutoff = cutoff
663 :
664 6550 : IF (multigrid_set) THEN
665 : ! Read the values from input
666 4 : IF (qs_control%commensurate_mgrids) THEN
667 0 : CPABORT("Do not specify cutoffs for the commensurate grids (NYI)")
668 : END IF
669 :
670 4 : CALL section_vals_val_get(mgrid_section, "MULTIGRID_CUTOFF", r_vals=cutofflist)
671 4 : IF (ASSOCIATED(cutofflist)) THEN
672 4 : IF (SIZE(cutofflist, 1) /= ngrid_level) THEN
673 0 : CPABORT("Number of multi-grids requested and number of cutoff values do not match")
674 : END IF
675 20 : DO igrid_level = 1, ngrid_level
676 20 : qs_control%e_cutoff(igrid_level) = cutofflist(igrid_level)
677 : END DO
678 : END IF
679 : ! set cutoff to smallest value in multgrid available with >= cutoff
680 20 : DO igrid_level = ngrid_level, 1, -1
681 16 : IF (qs_control%cutoff <= qs_control%e_cutoff(igrid_level)) THEN
682 0 : qs_control%cutoff = qs_control%e_cutoff(igrid_level)
683 0 : EXIT
684 : END IF
685 : ! set largest grid value to cutoff
686 20 : IF (igrid_level == 1) THEN
687 4 : qs_control%cutoff = qs_control%e_cutoff(1)
688 : END IF
689 : END DO
690 : ELSE
691 6546 : IF (qs_control%commensurate_mgrids) qs_control%progression_factor = 4.0_dp
692 6546 : qs_control%e_cutoff(1) = qs_control%cutoff
693 21698 : DO igrid_level = 2, ngrid_level
694 : qs_control%e_cutoff(igrid_level) = qs_control%e_cutoff(igrid_level - 1)/ &
695 21698 : qs_control%progression_factor
696 : END DO
697 : END IF
698 : ! check that multigrids are ordered
699 21714 : DO igrid_level = 2, ngrid_level
700 21714 : IF (qs_control%e_cutoff(igrid_level) > qs_control%e_cutoff(igrid_level - 1)) THEN
701 0 : CPABORT("The cutoff values for the multi-grids are not ordered from large to small")
702 15164 : ELSE IF (qs_control%e_cutoff(igrid_level) == qs_control%e_cutoff(igrid_level - 1)) THEN
703 0 : CPABORT("The same cutoff value was specified for two multi-grids")
704 : END IF
705 : END DO
706 6550 : CALL timestop(handle)
707 13100 : END SUBROUTINE read_mgrid_section
708 :
709 : ! **************************************************************************************************
710 : !> \brief ...
711 : !> \param qs_control ...
712 : !> \param qs_section ...
713 : ! **************************************************************************************************
714 104800 : SUBROUTINE read_qs_section(qs_control, qs_section)
715 :
716 : TYPE(qs_control_type), INTENT(INOUT) :: qs_control
717 : TYPE(section_vals_type), POINTER :: qs_section
718 :
719 : CHARACTER(len=*), PARAMETER :: routineN = 'read_qs_section'
720 :
721 : CHARACTER(LEN=default_string_length), &
722 6550 : DIMENSION(:), POINTER :: clist
723 : INTEGER :: handle, itmp, j, jj, k, n_rep, n_var, &
724 : ngauss, ngp, nrep
725 6550 : INTEGER, DIMENSION(:), POINTER :: tmplist
726 : LOGICAL :: explicit, was_present
727 : REAL(dp) :: tmp, tmpsqrt, value
728 6550 : REAL(dp), POINTER :: scal(:)
729 : TYPE(section_vals_type), POINTER :: cdft_control_section, ddapc_restraint_section, &
730 : dftb_parameter, dftb_section, genpot_section, lri_optbas_section, mull_section, &
731 : nonbonded_section, s2_restraint_section, se_section, xtb_parameter, xtb_section
732 :
733 6550 : CALL timeset(routineN, handle)
734 :
735 6550 : was_present = .FALSE.
736 6550 : NULLIFY (mull_section, ddapc_restraint_section, s2_restraint_section, &
737 6550 : se_section, dftb_section, xtb_section, dftb_parameter, xtb_parameter, lri_optbas_section, &
738 6550 : cdft_control_section, genpot_section)
739 :
740 6550 : mull_section => section_vals_get_subs_vals(qs_section, "MULLIKEN_RESTRAINT")
741 6550 : ddapc_restraint_section => section_vals_get_subs_vals(qs_section, "DDAPC_RESTRAINT")
742 6550 : s2_restraint_section => section_vals_get_subs_vals(qs_section, "S2_RESTRAINT")
743 6550 : se_section => section_vals_get_subs_vals(qs_section, "SE")
744 6550 : dftb_section => section_vals_get_subs_vals(qs_section, "DFTB")
745 6550 : xtb_section => section_vals_get_subs_vals(qs_section, "xTB")
746 6550 : dftb_parameter => section_vals_get_subs_vals(dftb_section, "PARAMETER")
747 6550 : xtb_parameter => section_vals_get_subs_vals(xtb_section, "PARAMETER")
748 6550 : lri_optbas_section => section_vals_get_subs_vals(qs_section, "OPTIMIZE_LRI_BASIS")
749 6550 : cdft_control_section => section_vals_get_subs_vals(qs_section, "CDFT")
750 6550 : nonbonded_section => section_vals_get_subs_vals(xtb_section, "NONBONDED")
751 6550 : genpot_section => section_vals_get_subs_vals(nonbonded_section, "GENPOT")
752 :
753 : ! Setup all defaults values and overwrite input parameters
754 : ! EPS_DEFAULT should set the target accuracy in the total energy (~per electron) or a closely related value
755 6550 : CALL section_vals_val_get(qs_section, "EPS_DEFAULT", r_val=value)
756 6550 : tmpsqrt = SQRT(value) ! a trick to work around a NAG 5.1 optimizer bug
757 :
758 : ! random choice ?
759 6550 : qs_control%eps_core_charge = value/100.0_dp
760 : ! correct if all Gaussians would have the same radius (overlap will be smaller than eps_pgf_orb**2).
761 : ! Can be significantly in error if not... requires fully new screening/pairlist procedures
762 6550 : qs_control%eps_pgf_orb = tmpsqrt
763 6550 : qs_control%eps_kg_orb = qs_control%eps_pgf_orb
764 : ! consistent since also a kind of overlap
765 6550 : qs_control%eps_ppnl = qs_control%eps_pgf_orb/100.0_dp
766 : ! accuracy is basically set by the overlap, this sets an empirical shift
767 6550 : qs_control%eps_ppl = 1.0E-2_dp
768 : !
769 6550 : qs_control%gapw_control%eps_cpc = value
770 : ! expexted error in the density
771 6550 : qs_control%eps_rho_gspace = value
772 6550 : qs_control%eps_rho_rspace = value
773 : ! error in the gradient, can be the sqrt of the error in the energy, ignored if map_consistent
774 6550 : qs_control%eps_gvg_rspace = tmpsqrt
775 : !
776 6550 : CALL section_vals_val_get(qs_section, "EPS_CORE_CHARGE", n_rep_val=n_rep)
777 6550 : IF (n_rep /= 0) THEN
778 0 : CALL section_vals_val_get(qs_section, "EPS_CORE_CHARGE", r_val=qs_control%eps_core_charge)
779 : END IF
780 6550 : CALL section_vals_val_get(qs_section, "EPS_GVG_RSPACE", n_rep_val=n_rep)
781 6550 : IF (n_rep /= 0) THEN
782 132 : CALL section_vals_val_get(qs_section, "EPS_GVG_RSPACE", r_val=qs_control%eps_gvg_rspace)
783 : END IF
784 6550 : CALL section_vals_val_get(qs_section, "EPS_PGF_ORB", n_rep_val=n_rep)
785 6550 : IF (n_rep /= 0) THEN
786 570 : CALL section_vals_val_get(qs_section, "EPS_PGF_ORB", r_val=qs_control%eps_pgf_orb)
787 : END IF
788 6550 : CALL section_vals_val_get(qs_section, "EPS_KG_ORB", n_rep_val=n_rep)
789 6550 : IF (n_rep /= 0) THEN
790 62 : CALL section_vals_val_get(qs_section, "EPS_KG_ORB", r_val=tmp)
791 62 : qs_control%eps_kg_orb = SQRT(tmp)
792 : END IF
793 6550 : CALL section_vals_val_get(qs_section, "EPS_PPL", n_rep_val=n_rep)
794 6550 : IF (n_rep /= 0) THEN
795 6550 : CALL section_vals_val_get(qs_section, "EPS_PPL", r_val=qs_control%eps_ppl)
796 : END IF
797 6550 : CALL section_vals_val_get(qs_section, "EPS_PPNL", n_rep_val=n_rep)
798 6550 : IF (n_rep /= 0) THEN
799 0 : CALL section_vals_val_get(qs_section, "EPS_PPNL", r_val=qs_control%eps_ppnl)
800 : END IF
801 6550 : CALL section_vals_val_get(qs_section, "EPS_RHO", n_rep_val=n_rep)
802 6550 : IF (n_rep /= 0) THEN
803 30 : CALL section_vals_val_get(qs_section, "EPS_RHO", r_val=qs_control%eps_rho_gspace)
804 30 : qs_control%eps_rho_rspace = qs_control%eps_rho_gspace
805 : END IF
806 6550 : CALL section_vals_val_get(qs_section, "EPS_RHO_RSPACE", n_rep_val=n_rep)
807 6550 : IF (n_rep /= 0) THEN
808 2 : CALL section_vals_val_get(qs_section, "EPS_RHO_RSPACE", r_val=qs_control%eps_rho_rspace)
809 : END IF
810 6550 : CALL section_vals_val_get(qs_section, "EPS_RHO_GSPACE", n_rep_val=n_rep)
811 6550 : IF (n_rep /= 0) THEN
812 2 : CALL section_vals_val_get(qs_section, "EPS_RHO_GSPACE", r_val=qs_control%eps_rho_gspace)
813 : END IF
814 6550 : CALL section_vals_val_get(qs_section, "EPS_FILTER_MATRIX", n_rep_val=n_rep)
815 6550 : IF (n_rep /= 0) THEN
816 6550 : CALL section_vals_val_get(qs_section, "EPS_FILTER_MATRIX", r_val=qs_control%eps_filter_matrix)
817 : END IF
818 6550 : CALL section_vals_val_get(qs_section, "EPS_CPC", n_rep_val=n_rep)
819 6550 : IF (n_rep /= 0) THEN
820 0 : CALL section_vals_val_get(qs_section, "EPS_CPC", r_val=qs_control%gapw_control%eps_cpc)
821 : END IF
822 :
823 6550 : CALL section_vals_val_get(qs_section, "EPSFIT", r_val=qs_control%gapw_control%eps_fit)
824 6550 : CALL section_vals_val_get(qs_section, "EPSISO", r_val=qs_control%gapw_control%eps_iso)
825 6550 : CALL section_vals_val_get(qs_section, "EPSSVD", r_val=qs_control%gapw_control%eps_svd)
826 6550 : CALL section_vals_val_get(qs_section, "EPSRHO0", r_val=qs_control%gapw_control%eps_Vrho0)
827 6550 : CALL section_vals_val_get(qs_section, "ALPHA0_HARD", r_val=qs_control%gapw_control%alpha0_hard)
828 6550 : qs_control%gapw_control%lrho1_eq_lrho0 = .FALSE.
829 6550 : qs_control%gapw_control%alpha0_hard_from_input = .FALSE.
830 6550 : IF (qs_control%gapw_control%alpha0_hard /= 0.0_dp) qs_control%gapw_control%alpha0_hard_from_input = .TRUE.
831 6550 : CALL section_vals_val_get(qs_section, "FORCE_PAW", l_val=qs_control%gapw_control%force_paw)
832 6550 : CALL section_vals_val_get(qs_section, "MAX_RAD_LOCAL", r_val=qs_control%gapw_control%max_rad_local)
833 :
834 6550 : CALL section_vals_val_get(qs_section, "MIN_PAIR_LIST_RADIUS", r_val=qs_control%pairlist_radius)
835 :
836 6550 : CALL section_vals_val_get(qs_section, "LS_SCF", l_val=qs_control%do_ls_scf)
837 6550 : CALL section_vals_val_get(qs_section, "ALMO_SCF", l_val=qs_control%do_almo_scf)
838 6550 : CALL section_vals_val_get(qs_section, "KG_METHOD", l_val=qs_control%do_kg)
839 :
840 : ! Logicals
841 6550 : CALL section_vals_val_get(qs_section, "REF_EMBED_SUBSYS", l_val=qs_control%ref_embed_subsys)
842 6550 : CALL section_vals_val_get(qs_section, "CLUSTER_EMBED_SUBSYS", l_val=qs_control%cluster_embed_subsys)
843 6550 : CALL section_vals_val_get(qs_section, "HIGH_LEVEL_EMBED_SUBSYS", l_val=qs_control%high_level_embed_subsys)
844 6550 : CALL section_vals_val_get(qs_section, "DFET_EMBEDDED", l_val=qs_control%dfet_embedded)
845 6550 : CALL section_vals_val_get(qs_section, "DMFET_EMBEDDED", l_val=qs_control%dmfet_embedded)
846 :
847 : ! Integers gapw
848 6550 : CALL section_vals_val_get(qs_section, "LMAXN1", i_val=qs_control%gapw_control%lmax_sphere)
849 6550 : CALL section_vals_val_get(qs_section, "LMAXN0", i_val=qs_control%gapw_control%lmax_rho0)
850 6550 : CALL section_vals_val_get(qs_section, "LADDN0", i_val=qs_control%gapw_control%ladd_rho0)
851 6550 : CALL section_vals_val_get(qs_section, "QUADRATURE", i_val=qs_control%gapw_control%quadrature)
852 : ! GAPW 1c basis
853 6550 : CALL section_vals_val_get(qs_section, "GAPW_1C_BASIS", i_val=qs_control%gapw_control%basis_1c)
854 6550 : IF (qs_control%gapw_control%basis_1c /= gapw_1c_orb) THEN
855 18 : qs_control%gapw_control%eps_svd = MAX(qs_control%gapw_control%eps_svd, 1.E-12_dp)
856 : END IF
857 :
858 : ! Integers grids
859 6550 : CALL section_vals_val_get(qs_section, "PW_GRID", i_val=itmp)
860 0 : SELECT CASE (itmp)
861 : CASE (do_pwgrid_spherical)
862 0 : qs_control%pw_grid_opt%spherical = .TRUE.
863 0 : qs_control%pw_grid_opt%fullspace = .FALSE.
864 : CASE (do_pwgrid_ns_fullspace)
865 6550 : qs_control%pw_grid_opt%spherical = .FALSE.
866 6550 : qs_control%pw_grid_opt%fullspace = .TRUE.
867 : CASE (do_pwgrid_ns_halfspace)
868 0 : qs_control%pw_grid_opt%spherical = .FALSE.
869 6550 : qs_control%pw_grid_opt%fullspace = .FALSE.
870 : END SELECT
871 :
872 : ! Method for PPL calculation
873 6550 : CALL section_vals_val_get(qs_section, "CORE_PPL", i_val=itmp)
874 6550 : qs_control%do_ppl_method = itmp
875 :
876 6550 : CALL section_vals_val_get(qs_section, "PW_GRID_LAYOUT", i_vals=tmplist)
877 19650 : qs_control%pw_grid_opt%distribution_layout = tmplist
878 6550 : CALL section_vals_val_get(qs_section, "PW_GRID_BLOCKED", i_val=qs_control%pw_grid_opt%blocked)
879 :
880 : !Integers extrapolation
881 6550 : CALL section_vals_val_get(qs_section, "EXTRAPOLATION", i_val=qs_control%wf_interpolation_method_nr)
882 6550 : CALL section_vals_val_get(qs_section, "EXTRAPOLATION_ORDER", i_val=qs_control%wf_extrapolation_order)
883 :
884 : !Method
885 6550 : CALL section_vals_val_get(qs_section, "METHOD", i_val=qs_control%method_id)
886 6550 : qs_control%gapw = .FALSE.
887 6550 : qs_control%gapw_xc = .FALSE.
888 6550 : qs_control%gpw = .FALSE.
889 6550 : qs_control%pao = .FALSE.
890 6550 : qs_control%dftb = .FALSE.
891 6550 : qs_control%xtb = .FALSE.
892 6550 : qs_control%semi_empirical = .FALSE.
893 6550 : qs_control%ofgpw = .FALSE.
894 6550 : qs_control%lrigpw = .FALSE.
895 6550 : qs_control%rigpw = .FALSE.
896 7356 : SELECT CASE (qs_control%method_id)
897 : CASE (do_method_gapw)
898 806 : CALL cite_reference(Lippert1999)
899 806 : CALL cite_reference(Krack2000)
900 806 : qs_control%gapw = .TRUE.
901 : CASE (do_method_gapw_xc)
902 106 : qs_control%gapw_xc = .TRUE.
903 : CASE (do_method_gpw)
904 4130 : CALL cite_reference(Lippert1997)
905 4130 : CALL cite_reference(VandeVondele2005a)
906 4130 : qs_control%gpw = .TRUE.
907 : CASE (do_method_ofgpw)
908 0 : qs_control%ofgpw = .TRUE.
909 : CASE (do_method_lrigpw)
910 40 : qs_control%lrigpw = .TRUE.
911 : CASE (do_method_rigpw)
912 0 : qs_control%rigpw = .TRUE.
913 : CASE (do_method_dftb)
914 222 : qs_control%dftb = .TRUE.
915 222 : CALL cite_reference(Porezag1995)
916 222 : CALL cite_reference(Seifert1996)
917 : CASE (do_method_xtb)
918 248 : qs_control%xtb = .TRUE.
919 248 : CALL cite_reference(Grimme2017)
920 : CASE (do_method_mndo)
921 52 : CALL cite_reference(Dewar1977)
922 52 : qs_control%semi_empirical = .TRUE.
923 : CASE (do_method_am1)
924 112 : CALL cite_reference(Dewar1985)
925 112 : qs_control%semi_empirical = .TRUE.
926 : CASE (do_method_pm3)
927 46 : CALL cite_reference(Stewart1989)
928 46 : qs_control%semi_empirical = .TRUE.
929 : CASE (do_method_pnnl)
930 14 : CALL cite_reference(Schenter2008)
931 14 : qs_control%semi_empirical = .TRUE.
932 : CASE (do_method_pm6)
933 754 : CALL cite_reference(Stewart2007)
934 754 : qs_control%semi_empirical = .TRUE.
935 : CASE (do_method_pm6fm)
936 0 : CALL cite_reference(VanVoorhis2015)
937 0 : qs_control%semi_empirical = .TRUE.
938 : CASE (do_method_pdg)
939 2 : CALL cite_reference(Repasky2002)
940 2 : qs_control%semi_empirical = .TRUE.
941 : CASE (do_method_rm1)
942 2 : CALL cite_reference(Rocha2006)
943 2 : qs_control%semi_empirical = .TRUE.
944 : CASE (do_method_mndod)
945 16 : CALL cite_reference(Dewar1977)
946 16 : CALL cite_reference(Thiel1992)
947 6566 : qs_control%semi_empirical = .TRUE.
948 : END SELECT
949 :
950 6550 : CALL section_vals_get(mull_section, explicit=qs_control%mulliken_restraint)
951 :
952 6550 : IF (qs_control%mulliken_restraint) THEN
953 2 : CALL section_vals_val_get(mull_section, "STRENGTH", r_val=qs_control%mulliken_restraint_control%strength)
954 2 : CALL section_vals_val_get(mull_section, "TARGET", r_val=qs_control%mulliken_restraint_control%target)
955 2 : CALL section_vals_val_get(mull_section, "ATOMS", n_rep_val=n_rep)
956 2 : jj = 0
957 4 : DO k = 1, n_rep
958 2 : CALL section_vals_val_get(mull_section, "ATOMS", i_rep_val=k, i_vals=tmplist)
959 4 : jj = jj + SIZE(tmplist)
960 : END DO
961 2 : qs_control%mulliken_restraint_control%natoms = jj
962 2 : IF (qs_control%mulliken_restraint_control%natoms < 1) &
963 0 : CPABORT("Need at least 1 atom to use mulliken constraints")
964 6 : ALLOCATE (qs_control%mulliken_restraint_control%atoms(qs_control%mulliken_restraint_control%natoms))
965 2 : jj = 0
966 6 : DO k = 1, n_rep
967 2 : CALL section_vals_val_get(mull_section, "ATOMS", i_rep_val=k, i_vals=tmplist)
968 6 : DO j = 1, SIZE(tmplist)
969 2 : jj = jj + 1
970 4 : qs_control%mulliken_restraint_control%atoms(jj) = tmplist(j)
971 : END DO
972 : END DO
973 : END IF
974 6550 : CALL section_vals_get(ddapc_restraint_section, n_repetition=nrep, explicit=qs_control%ddapc_restraint)
975 6550 : IF (qs_control%ddapc_restraint) THEN
976 60 : ALLOCATE (qs_control%ddapc_restraint_control(nrep))
977 14 : CALL read_ddapc_section(qs_control, qs_section=qs_section)
978 14 : qs_control%ddapc_restraint_is_spin = .FALSE.
979 14 : qs_control%ddapc_explicit_potential = .FALSE.
980 : END IF
981 :
982 6550 : CALL section_vals_get(s2_restraint_section, explicit=qs_control%s2_restraint)
983 6550 : IF (qs_control%s2_restraint) THEN
984 : CALL section_vals_val_get(s2_restraint_section, "STRENGTH", &
985 0 : r_val=qs_control%s2_restraint_control%strength)
986 : CALL section_vals_val_get(s2_restraint_section, "TARGET", &
987 0 : r_val=qs_control%s2_restraint_control%target)
988 : CALL section_vals_val_get(s2_restraint_section, "FUNCTIONAL_FORM", &
989 0 : i_val=qs_control%s2_restraint_control%functional_form)
990 : END IF
991 :
992 6550 : CALL section_vals_get(cdft_control_section, explicit=qs_control%cdft)
993 6550 : IF (qs_control%cdft) THEN
994 264 : CALL read_cdft_control_section(qs_control, cdft_control_section)
995 : END IF
996 :
997 : ! Semi-empirical code
998 6550 : IF (qs_control%semi_empirical) THEN
999 : CALL section_vals_val_get(se_section, "ORTHOGONAL_BASIS", &
1000 998 : l_val=qs_control%se_control%orthogonal_basis)
1001 : CALL section_vals_val_get(se_section, "DELTA", &
1002 998 : r_val=qs_control%se_control%delta)
1003 : CALL section_vals_val_get(se_section, "ANALYTICAL_GRADIENTS", &
1004 998 : l_val=qs_control%se_control%analytical_gradients)
1005 : CALL section_vals_val_get(se_section, "FORCE_KDSO-D_EXCHANGE", &
1006 998 : l_val=qs_control%se_control%force_kdsod_EX)
1007 : ! Integral Screening
1008 : CALL section_vals_val_get(se_section, "INTEGRAL_SCREENING", &
1009 998 : i_val=qs_control%se_control%integral_screening)
1010 998 : IF (qs_control%method_id == do_method_pnnl) THEN
1011 14 : IF (qs_control%se_control%integral_screening /= do_se_IS_slater) &
1012 : CALL cp_warn(__LOCATION__, &
1013 : "PNNL semi-empirical parameterization supports only the Slater type "// &
1014 0 : "integral scheme. Revert to Slater and continue the calculation.")
1015 14 : qs_control%se_control%integral_screening = do_se_IS_slater
1016 : END IF
1017 : ! Global Arrays variable
1018 : CALL section_vals_val_get(se_section, "GA%NCELLS", &
1019 998 : i_val=qs_control%se_control%ga_ncells)
1020 : ! Long-Range correction
1021 : CALL section_vals_val_get(se_section, "LR_CORRECTION%CUTOFF", &
1022 998 : r_val=qs_control%se_control%cutoff_lrc)
1023 998 : qs_control%se_control%taper_lrc = qs_control%se_control%cutoff_lrc
1024 : CALL section_vals_val_get(se_section, "LR_CORRECTION%RC_TAPER", &
1025 998 : explicit=explicit)
1026 998 : IF (explicit) THEN
1027 : CALL section_vals_val_get(se_section, "LR_CORRECTION%RC_TAPER", &
1028 0 : r_val=qs_control%se_control%taper_lrc)
1029 : END IF
1030 : CALL section_vals_val_get(se_section, "LR_CORRECTION%RC_RANGE", &
1031 998 : r_val=qs_control%se_control%range_lrc)
1032 : ! Coulomb
1033 : CALL section_vals_val_get(se_section, "COULOMB%CUTOFF", &
1034 998 : r_val=qs_control%se_control%cutoff_cou)
1035 998 : qs_control%se_control%taper_cou = qs_control%se_control%cutoff_cou
1036 : CALL section_vals_val_get(se_section, "COULOMB%RC_TAPER", &
1037 998 : explicit=explicit)
1038 998 : IF (explicit) THEN
1039 : CALL section_vals_val_get(se_section, "COULOMB%RC_TAPER", &
1040 0 : r_val=qs_control%se_control%taper_cou)
1041 : END IF
1042 : CALL section_vals_val_get(se_section, "COULOMB%RC_RANGE", &
1043 998 : r_val=qs_control%se_control%range_cou)
1044 : ! Exchange
1045 : CALL section_vals_val_get(se_section, "EXCHANGE%CUTOFF", &
1046 998 : r_val=qs_control%se_control%cutoff_exc)
1047 998 : qs_control%se_control%taper_exc = qs_control%se_control%cutoff_exc
1048 : CALL section_vals_val_get(se_section, "EXCHANGE%RC_TAPER", &
1049 998 : explicit=explicit)
1050 998 : IF (explicit) THEN
1051 : CALL section_vals_val_get(se_section, "EXCHANGE%RC_TAPER", &
1052 38 : r_val=qs_control%se_control%taper_exc)
1053 : END IF
1054 : CALL section_vals_val_get(se_section, "EXCHANGE%RC_RANGE", &
1055 998 : r_val=qs_control%se_control%range_exc)
1056 : ! Screening (only if the integral scheme is of dumped type)
1057 998 : IF (qs_control%se_control%integral_screening == do_se_IS_kdso_d) THEN
1058 : CALL section_vals_val_get(se_section, "SCREENING%RC_TAPER", &
1059 14 : r_val=qs_control%se_control%taper_scr)
1060 : CALL section_vals_val_get(se_section, "SCREENING%RC_RANGE", &
1061 14 : r_val=qs_control%se_control%range_scr)
1062 : END IF
1063 : ! Periodic Type Calculation
1064 : CALL section_vals_val_get(se_section, "PERIODIC", &
1065 998 : i_val=qs_control%se_control%periodic_type)
1066 1964 : SELECT CASE (qs_control%se_control%periodic_type)
1067 : CASE (do_se_lr_none)
1068 966 : qs_control%se_control%do_ewald = .FALSE.
1069 966 : qs_control%se_control%do_ewald_r3 = .FALSE.
1070 966 : qs_control%se_control%do_ewald_gks = .FALSE.
1071 : CASE (do_se_lr_ewald)
1072 30 : qs_control%se_control%do_ewald = .TRUE.
1073 30 : qs_control%se_control%do_ewald_r3 = .FALSE.
1074 30 : qs_control%se_control%do_ewald_gks = .FALSE.
1075 : CASE (do_se_lr_ewald_gks)
1076 2 : qs_control%se_control%do_ewald = .FALSE.
1077 2 : qs_control%se_control%do_ewald_r3 = .FALSE.
1078 2 : qs_control%se_control%do_ewald_gks = .TRUE.
1079 2 : IF (qs_control%method_id /= do_method_pnnl) &
1080 : CALL cp_abort(__LOCATION__, &
1081 : "A periodic semi-empirical calculation was requested with a long-range "// &
1082 : "summation on the single integral evaluation. This scheme is supported "// &
1083 0 : "only by the PNNL parameterization.")
1084 : CASE (do_se_lr_ewald_r3)
1085 0 : qs_control%se_control%do_ewald = .TRUE.
1086 0 : qs_control%se_control%do_ewald_r3 = .TRUE.
1087 0 : qs_control%se_control%do_ewald_gks = .FALSE.
1088 0 : IF (qs_control%se_control%integral_screening /= do_se_IS_kdso) &
1089 : CALL cp_abort(__LOCATION__, &
1090 : "A periodic semi-empirical calculation was requested with a long-range "// &
1091 : "summation for the slowly convergent part 1/R^3, which is not congruent "// &
1092 : "with the integral screening chosen. The only integral screening supported "// &
1093 998 : "by this periodic type calculation is the standard Klopman-Dewar-Sabelli-Ohno.")
1094 : END SELECT
1095 :
1096 : ! dispersion pair potentials
1097 : CALL section_vals_val_get(se_section, "DISPERSION", &
1098 998 : l_val=qs_control%se_control%dispersion)
1099 : CALL section_vals_val_get(se_section, "DISPERSION_RADIUS", &
1100 998 : r_val=qs_control%se_control%rcdisp)
1101 : CALL section_vals_val_get(se_section, "COORDINATION_CUTOFF", &
1102 998 : r_val=qs_control%se_control%epscn)
1103 998 : CALL section_vals_val_get(se_section, "D3_SCALING", r_vals=scal)
1104 998 : qs_control%se_control%sd3(1) = scal(1)
1105 998 : qs_control%se_control%sd3(2) = scal(2)
1106 998 : qs_control%se_control%sd3(3) = scal(3)
1107 : CALL section_vals_val_get(se_section, "DISPERSION_PARAMETER_FILE", &
1108 998 : c_val=qs_control%se_control%dispersion_parameter_file)
1109 :
1110 : ! Stop the execution for non-implemented features
1111 998 : IF (qs_control%se_control%periodic_type == do_se_lr_ewald_r3) THEN
1112 0 : CPABORT("EWALD_R3 not implemented yet!")
1113 : END IF
1114 :
1115 : IF (qs_control%method_id == do_method_mndo .OR. &
1116 : qs_control%method_id == do_method_am1 .OR. &
1117 : qs_control%method_id == do_method_mndod .OR. &
1118 : qs_control%method_id == do_method_pdg .OR. &
1119 : qs_control%method_id == do_method_pm3 .OR. &
1120 : qs_control%method_id == do_method_pm6 .OR. &
1121 : qs_control%method_id == do_method_pm6fm .OR. &
1122 998 : qs_control%method_id == do_method_pnnl .OR. &
1123 : qs_control%method_id == do_method_rm1) THEN
1124 998 : qs_control%se_control%orthogonal_basis = .TRUE.
1125 : END IF
1126 : END IF
1127 :
1128 : ! DFTB code
1129 6550 : IF (qs_control%dftb) THEN
1130 : CALL section_vals_val_get(dftb_section, "ORTHOGONAL_BASIS", &
1131 222 : l_val=qs_control%dftb_control%orthogonal_basis)
1132 : CALL section_vals_val_get(dftb_section, "SELF_CONSISTENT", &
1133 222 : l_val=qs_control%dftb_control%self_consistent)
1134 : CALL section_vals_val_get(dftb_section, "DISPERSION", &
1135 222 : l_val=qs_control%dftb_control%dispersion)
1136 : CALL section_vals_val_get(dftb_section, "DIAGONAL_DFTB3", &
1137 222 : l_val=qs_control%dftb_control%dftb3_diagonal)
1138 : CALL section_vals_val_get(dftb_section, "HB_SR_GAMMA", &
1139 222 : l_val=qs_control%dftb_control%hb_sr_damp)
1140 : CALL section_vals_val_get(dftb_section, "EPS_DISP", &
1141 222 : r_val=qs_control%dftb_control%eps_disp)
1142 222 : CALL section_vals_val_get(dftb_section, "DO_EWALD", explicit=explicit)
1143 222 : IF (explicit) THEN
1144 : CALL section_vals_val_get(dftb_section, "DO_EWALD", &
1145 166 : l_val=qs_control%dftb_control%do_ewald)
1146 : ELSE
1147 56 : qs_control%dftb_control%do_ewald = (qs_control%periodicity /= 0)
1148 : END IF
1149 : CALL section_vals_val_get(dftb_parameter, "PARAM_FILE_PATH", &
1150 222 : c_val=qs_control%dftb_control%sk_file_path)
1151 : CALL section_vals_val_get(dftb_parameter, "PARAM_FILE_NAME", &
1152 222 : c_val=qs_control%dftb_control%sk_file_list)
1153 : CALL section_vals_val_get(dftb_parameter, "HB_SR_PARAM", &
1154 222 : r_val=qs_control%dftb_control%hb_sr_para)
1155 222 : CALL section_vals_val_get(dftb_parameter, "SK_FILE", n_rep_val=n_var)
1156 470 : ALLOCATE (qs_control%dftb_control%sk_pair_list(3, n_var))
1157 284 : DO k = 1, n_var
1158 : CALL section_vals_val_get(dftb_parameter, "SK_FILE", i_rep_val=k, &
1159 62 : c_vals=clist)
1160 470 : qs_control%dftb_control%sk_pair_list(1:3, k) = clist(1:3)
1161 : END DO
1162 : ! Dispersion type
1163 : CALL section_vals_val_get(dftb_parameter, "DISPERSION_TYPE", &
1164 222 : i_val=qs_control%dftb_control%dispersion_type)
1165 : CALL section_vals_val_get(dftb_parameter, "UFF_FORCE_FIELD", &
1166 222 : c_val=qs_control%dftb_control%uff_force_field)
1167 : ! D3 Dispersion
1168 : CALL section_vals_val_get(dftb_parameter, "DISPERSION_RADIUS", &
1169 222 : r_val=qs_control%dftb_control%rcdisp)
1170 : CALL section_vals_val_get(dftb_parameter, "COORDINATION_CUTOFF", &
1171 222 : r_val=qs_control%dftb_control%epscn)
1172 : CALL section_vals_val_get(dftb_parameter, "D2_EXP_PRE", &
1173 222 : r_val=qs_control%dftb_control%exp_pre)
1174 : CALL section_vals_val_get(dftb_parameter, "D2_SCALING", &
1175 222 : r_val=qs_control%dftb_control%scaling)
1176 222 : CALL section_vals_val_get(dftb_parameter, "D3_SCALING", r_vals=scal)
1177 222 : qs_control%dftb_control%sd3(1) = scal(1)
1178 222 : qs_control%dftb_control%sd3(2) = scal(2)
1179 222 : qs_control%dftb_control%sd3(3) = scal(3)
1180 222 : CALL section_vals_val_get(dftb_parameter, "D3BJ_SCALING", r_vals=scal)
1181 222 : qs_control%dftb_control%sd3bj(1) = scal(1)
1182 222 : qs_control%dftb_control%sd3bj(2) = scal(2)
1183 222 : qs_control%dftb_control%sd3bj(3) = scal(3)
1184 222 : qs_control%dftb_control%sd3bj(4) = scal(4)
1185 : CALL section_vals_val_get(dftb_parameter, "DISPERSION_PARAMETER_FILE", &
1186 222 : c_val=qs_control%dftb_control%dispersion_parameter_file)
1187 :
1188 222 : IF (qs_control%dftb_control%dispersion) CALL cite_reference(Zhechkov2005)
1189 222 : IF (qs_control%dftb_control%self_consistent) CALL cite_reference(Elstner1998)
1190 666 : IF (qs_control%dftb_control%hb_sr_damp) CALL cite_reference(Hu2007)
1191 : END IF
1192 :
1193 : ! xTB code
1194 6550 : IF (qs_control%xtb) THEN
1195 248 : CALL section_vals_val_get(xtb_section, "DO_EWALD", explicit=explicit)
1196 248 : IF (explicit) THEN
1197 : CALL section_vals_val_get(xtb_section, "DO_EWALD", &
1198 134 : l_val=qs_control%xtb_control%do_ewald)
1199 : ELSE
1200 114 : qs_control%xtb_control%do_ewald = (qs_control%periodicity /= 0)
1201 : END IF
1202 248 : CALL section_vals_val_get(xtb_section, "STO_NG", i_val=ngauss)
1203 248 : qs_control%xtb_control%sto_ng = ngauss
1204 248 : CALL section_vals_val_get(xtb_section, "HYDROGEN_STO_NG", i_val=ngauss)
1205 248 : qs_control%xtb_control%h_sto_ng = ngauss
1206 : CALL section_vals_val_get(xtb_parameter, "PARAM_FILE_PATH", &
1207 248 : c_val=qs_control%xtb_control%parameter_file_path)
1208 : CALL section_vals_val_get(xtb_parameter, "PARAM_FILE_NAME", &
1209 248 : c_val=qs_control%xtb_control%parameter_file_name)
1210 : ! D3 Dispersion
1211 : CALL section_vals_val_get(xtb_parameter, "DISPERSION_RADIUS", &
1212 248 : r_val=qs_control%xtb_control%rcdisp)
1213 : CALL section_vals_val_get(xtb_parameter, "COORDINATION_CUTOFF", &
1214 248 : r_val=qs_control%xtb_control%epscn)
1215 248 : CALL section_vals_val_get(xtb_parameter, "D3BJ_SCALING", r_vals=scal)
1216 248 : qs_control%xtb_control%s6 = scal(1)
1217 248 : qs_control%xtb_control%s8 = scal(2)
1218 248 : CALL section_vals_val_get(xtb_parameter, "D3BJ_PARAM", r_vals=scal)
1219 248 : qs_control%xtb_control%a1 = scal(1)
1220 248 : qs_control%xtb_control%a2 = scal(2)
1221 : CALL section_vals_val_get(xtb_parameter, "DISPERSION_PARAMETER_FILE", &
1222 248 : c_val=qs_control%xtb_control%dispersion_parameter_file)
1223 : ! global parameters
1224 248 : CALL section_vals_val_get(xtb_parameter, "HUCKEL_CONSTANTS", r_vals=scal)
1225 248 : qs_control%xtb_control%ks = scal(1)
1226 248 : qs_control%xtb_control%kp = scal(2)
1227 248 : qs_control%xtb_control%kd = scal(3)
1228 248 : qs_control%xtb_control%ksp = scal(4)
1229 248 : qs_control%xtb_control%k2sh = scal(5)
1230 248 : CALL section_vals_val_get(xtb_parameter, "COULOMB_CONSTANTS", r_vals=scal)
1231 248 : qs_control%xtb_control%kg = scal(1)
1232 248 : qs_control%xtb_control%kf = scal(2)
1233 248 : CALL section_vals_val_get(xtb_parameter, "CN_CONSTANTS", r_vals=scal)
1234 248 : qs_control%xtb_control%kcns = scal(1)
1235 248 : qs_control%xtb_control%kcnp = scal(2)
1236 248 : qs_control%xtb_control%kcnd = scal(3)
1237 248 : CALL section_vals_val_get(xtb_parameter, "EN_CONSTANT", r_vals=scal)
1238 248 : qs_control%xtb_control%ken = scal(1)
1239 : ! XB
1240 : CALL section_vals_val_get(xtb_section, "USE_HALOGEN_CORRECTION", &
1241 248 : l_val=qs_control%xtb_control%xb_interaction)
1242 248 : CALL section_vals_val_get(xtb_parameter, "HALOGEN_BINDING", r_vals=scal)
1243 248 : qs_control%xtb_control%kxr = scal(1)
1244 248 : qs_control%xtb_control%kx2 = scal(2)
1245 : ! NONBONDED interactions
1246 : CALL section_vals_val_get(xtb_section, "DO_NONBONDED", &
1247 248 : l_val=qs_control%xtb_control%do_nonbonded)
1248 248 : CALL section_vals_get(nonbonded_section, explicit=explicit)
1249 248 : IF (explicit .AND. qs_control%xtb_control%do_nonbonded) THEN
1250 6 : CALL section_vals_get(genpot_section, explicit=explicit, n_repetition=ngp)
1251 6 : IF (explicit) THEN
1252 6 : CALL pair_potential_reallocate(qs_control%xtb_control%nonbonded, 1, ngp, gp=.TRUE.)
1253 6 : CALL read_gp_section(qs_control%xtb_control%nonbonded, genpot_section, 0)
1254 : END IF
1255 : END IF !nonbonded
1256 : ! SR Coulomb
1257 : CALL section_vals_val_get(xtb_section, "OLD_COULOMB_DAMPING", &
1258 248 : l_val=qs_control%xtb_control%old_coulomb_damping)
1259 248 : CALL section_vals_val_get(xtb_parameter, "COULOMB_SR_CUT", r_vals=scal)
1260 248 : qs_control%xtb_control%coulomb_sr_cut = scal(1)
1261 248 : CALL section_vals_val_get(xtb_parameter, "COULOMB_SR_EPS", r_vals=scal)
1262 248 : qs_control%xtb_control%coulomb_sr_eps = scal(1)
1263 : ! XB_radius
1264 248 : CALL section_vals_val_get(xtb_parameter, "XB_RADIUS", r_val=qs_control%xtb_control%xb_radius)
1265 : ! Kab
1266 248 : CALL section_vals_val_get(xtb_parameter, "KAB_PARAM", n_rep_val=n_rep)
1267 : ! For debug purposes
1268 : CALL section_vals_val_get(xtb_section, "COULOMB_INTERACTION", &
1269 248 : l_val=qs_control%xtb_control%coulomb_interaction)
1270 : CALL section_vals_val_get(xtb_section, "COULOMB_LR", &
1271 248 : l_val=qs_control%xtb_control%coulomb_lr)
1272 : CALL section_vals_val_get(xtb_section, "TB3_INTERACTION", &
1273 248 : l_val=qs_control%xtb_control%tb3_interaction)
1274 : ! Check for bad atomic charges
1275 : CALL section_vals_val_get(xtb_section, "CHECK_ATOMIC_CHARGES", &
1276 248 : l_val=qs_control%xtb_control%check_atomic_charges)
1277 :
1278 248 : qs_control%xtb_control%kab_nval = n_rep
1279 248 : IF (n_rep > 0) THEN
1280 6 : ALLOCATE (qs_control%xtb_control%kab_param(3, n_rep))
1281 6 : ALLOCATE (qs_control%xtb_control%kab_types(2, n_rep))
1282 6 : ALLOCATE (qs_control%xtb_control%kab_vals(n_rep))
1283 4 : DO j = 1, n_rep
1284 2 : CALL section_vals_val_get(xtb_parameter, "KAB_PARAM", i_rep_val=j, c_vals=clist)
1285 2 : qs_control%xtb_control%kab_param(1, j) = clist(1)
1286 : CALL get_ptable_info(clist(1), &
1287 2 : ielement=qs_control%xtb_control%kab_types(1, j))
1288 2 : qs_control%xtb_control%kab_param(2, j) = clist(2)
1289 : CALL get_ptable_info(clist(2), &
1290 2 : ielement=qs_control%xtb_control%kab_types(2, j))
1291 2 : qs_control%xtb_control%kab_param(3, j) = clist(3)
1292 4 : READ (clist(3), '(F10.0)') qs_control%xtb_control%kab_vals(j)
1293 : END DO
1294 : END IF
1295 : END IF
1296 :
1297 : ! Optimize LRI basis set
1298 6550 : CALL section_vals_get(lri_optbas_section, explicit=qs_control%lri_optbas)
1299 :
1300 6550 : CALL timestop(handle)
1301 6550 : END SUBROUTINE read_qs_section
1302 :
1303 : ! **************************************************************************************************
1304 : !> \brief ...
1305 : !> \param t_control ...
1306 : !> \param dft_section ...
1307 : ! **************************************************************************************************
1308 12 : SUBROUTINE read_tddfpt_control(t_control, dft_section)
1309 : TYPE(tddfpt_control_type) :: t_control
1310 : TYPE(section_vals_type), POINTER :: dft_section
1311 :
1312 : LOGICAL :: kenergy_den
1313 : TYPE(section_vals_type), POINTER :: sic_section, t_section
1314 :
1315 12 : kenergy_den = .FALSE.
1316 12 : NULLIFY (sic_section, t_section)
1317 12 : t_section => section_vals_get_subs_vals(dft_section, "TDDFPT")
1318 :
1319 12 : CALL section_vals_val_get(t_section, "CONVERGENCE", r_val=t_control%tolerance)
1320 12 : CALL section_vals_val_get(t_section, "NEV", i_val=t_control%n_ev)
1321 12 : CALL section_vals_val_get(t_section, "MAX_KV", i_val=t_control%max_kv)
1322 12 : CALL section_vals_val_get(t_section, "RESTARTS", i_val=t_control%n_restarts)
1323 12 : CALL section_vals_val_get(t_section, "NREORTHO", i_val=t_control%n_reortho)
1324 12 : CALL section_vals_val_get(t_section, "RES_ETYPE", i_val=t_control%res_etype)
1325 12 : CALL section_vals_val_get(t_section, "DIAG_METHOD", i_val=t_control%diag_method)
1326 12 : CALL section_vals_val_get(t_section, "KERNEL", l_val=t_control%do_kernel)
1327 12 : CALL section_vals_val_get(t_section, "LSD_SINGLETS", l_val=t_control%lsd_singlets)
1328 12 : CALL section_vals_val_get(t_section, "INVERT_S", l_val=t_control%invert_S)
1329 12 : CALL section_vals_val_get(t_section, "PRECOND", l_val=t_control%precond)
1330 12 : CALL section_vals_val_get(t_section, "OE_CORR", i_val=t_control%oe_corr)
1331 :
1332 12 : t_control%use_kinetic_energy_density = .FALSE.
1333 12 : sic_section => section_vals_get_subs_vals(t_section, "SIC")
1334 12 : CALL section_vals_val_get(sic_section, "SIC_METHOD", i_val=t_control%sic_method_id)
1335 12 : CALL section_vals_val_get(sic_section, "ORBITAL_SET", i_val=t_control%sic_list_id)
1336 12 : CALL section_vals_val_get(sic_section, "SIC_SCALING_A", r_val=t_control%sic_scaling_a)
1337 12 : CALL section_vals_val_get(sic_section, "SIC_SCALING_B", r_val=t_control%sic_scaling_b)
1338 :
1339 12 : END SUBROUTINE read_tddfpt_control
1340 :
1341 : ! **************************************************************************************************
1342 : !> \brief Read TDDFPT-related input parameters.
1343 : !> \param t_control TDDFPT control parameters
1344 : !> \param t_section TDDFPT input section
1345 : !> \param qs_control Quickstep control parameters
1346 : ! **************************************************************************************************
1347 6550 : SUBROUTINE read_tddfpt2_control(t_control, t_section, qs_control)
1348 : TYPE(tddfpt2_control_type), POINTER :: t_control
1349 : TYPE(section_vals_type), POINTER :: t_section
1350 : TYPE(qs_control_type), POINTER :: qs_control
1351 :
1352 : CHARACTER(LEN=*), PARAMETER :: routineN = 'read_tddfpt2_control'
1353 :
1354 : CHARACTER(LEN=default_string_length), &
1355 6550 : DIMENSION(:), POINTER :: tmpstringlist
1356 : INTEGER :: handle, irep, isize, nrep
1357 6550 : INTEGER, ALLOCATABLE, DIMENSION(:) :: inds
1358 : LOGICAL :: do_ewald, do_exchange, expl, explicit, &
1359 : multigrid_set
1360 : REAL(KIND=dp) :: filter, fval, hfx
1361 : TYPE(section_vals_type), POINTER :: dipole_section, mgrid_section, &
1362 : soc_section, stda_section, xc_func, &
1363 : xc_section
1364 :
1365 6550 : CALL timeset(routineN, handle)
1366 :
1367 6550 : CALL section_vals_val_get(t_section, "_SECTION_PARAMETERS_", l_val=t_control%enabled)
1368 :
1369 6550 : CALL section_vals_val_get(t_section, "NSTATES", i_val=t_control%nstates)
1370 6550 : CALL section_vals_val_get(t_section, "MAX_ITER", i_val=t_control%niters)
1371 6550 : CALL section_vals_val_get(t_section, "MAX_KV", i_val=t_control%nkvs)
1372 6550 : CALL section_vals_val_get(t_section, "NLUMO", i_val=t_control%nlumo)
1373 6550 : CALL section_vals_val_get(t_section, "NPROC_STATE", i_val=t_control%nprocs)
1374 6550 : CALL section_vals_val_get(t_section, "KERNEL", i_val=t_control%kernel)
1375 6550 : CALL section_vals_val_get(t_section, "OE_CORR", i_val=t_control%oe_corr)
1376 6550 : CALL section_vals_val_get(t_section, "EV_SHIFT", r_val=t_control%ev_shift)
1377 6550 : CALL section_vals_val_get(t_section, "EOS_SHIFT", r_val=t_control%eos_shift)
1378 :
1379 6550 : CALL section_vals_val_get(t_section, "CONVERGENCE", r_val=t_control%conv)
1380 6550 : CALL section_vals_val_get(t_section, "MIN_AMPLITUDE", r_val=t_control%min_excitation_amplitude)
1381 6550 : CALL section_vals_val_get(t_section, "ORTHOGONAL_EPS", r_val=t_control%orthogonal_eps)
1382 :
1383 6550 : CALL section_vals_val_get(t_section, "RESTART", l_val=t_control%is_restart)
1384 6550 : CALL section_vals_val_get(t_section, "RKS_TRIPLETS", l_val=t_control%rks_triplets)
1385 6550 : CALL section_vals_val_get(t_section, "DO_LRIGPW", l_val=t_control%do_lrigpw)
1386 6550 : CALL section_vals_val_get(t_section, "ADMM_KERNEL_CORRECTION_SYMMETRIC", l_val=t_control%admm_symm)
1387 6550 : CALL section_vals_val_get(t_section, "ADMM_KERNEL_XC_CORRECTION", l_val=t_control%admm_xc_correction)
1388 :
1389 : ! read automatically generated auxiliary basis for LRI
1390 6550 : CALL section_vals_val_get(t_section, "AUTO_BASIS", n_rep_val=nrep)
1391 13100 : DO irep = 1, nrep
1392 6550 : CALL section_vals_val_get(t_section, "AUTO_BASIS", i_rep_val=irep, c_vals=tmpstringlist)
1393 13100 : IF (SIZE(tmpstringlist) == 2) THEN
1394 6550 : CALL uppercase(tmpstringlist(2))
1395 6550 : SELECT CASE (tmpstringlist(2))
1396 : CASE ("X")
1397 0 : isize = -1
1398 : CASE ("SMALL")
1399 0 : isize = 0
1400 : CASE ("MEDIUM")
1401 0 : isize = 1
1402 : CASE ("LARGE")
1403 0 : isize = 2
1404 : CASE ("HUGE")
1405 0 : isize = 3
1406 : CASE DEFAULT
1407 6550 : CPABORT("Unknown basis size in AUTO_BASIS keyword:"//TRIM(tmpstringlist(1)))
1408 : END SELECT
1409 : !
1410 6550 : SELECT CASE (tmpstringlist(1))
1411 : CASE ("X")
1412 : CASE ("P_LRI_AUX")
1413 0 : t_control%auto_basis_p_lri_aux = isize
1414 : CASE DEFAULT
1415 6550 : CPABORT("Unknown basis type in AUTO_BASIS keyword:"//TRIM(tmpstringlist(1)))
1416 : END SELECT
1417 : ELSE
1418 : CALL cp_abort(__LOCATION__, &
1419 0 : "AUTO_BASIS keyword in &PROPERTIES &TDDFT section has a wrong number of arguments.")
1420 : END IF
1421 : END DO
1422 :
1423 6550 : IF (t_control%conv < 0) &
1424 0 : t_control%conv = ABS(t_control%conv)
1425 :
1426 : ! DIPOLE_MOMENTS subsection
1427 6550 : dipole_section => section_vals_get_subs_vals(t_section, "DIPOLE_MOMENTS")
1428 6550 : CALL section_vals_val_get(dipole_section, "DIPOLE_FORM", explicit=explicit)
1429 6550 : IF (explicit) THEN
1430 10 : CALL section_vals_val_get(dipole_section, "DIPOLE_FORM", i_val=t_control%dipole_form)
1431 : ELSE
1432 6540 : t_control%dipole_form = 0
1433 : END IF
1434 6550 : CALL section_vals_val_get(dipole_section, "REFERENCE", i_val=t_control%dipole_reference)
1435 6550 : CALL section_vals_val_get(dipole_section, "REFERENCE_POINT", explicit=explicit)
1436 6550 : IF (explicit) THEN
1437 0 : CALL section_vals_val_get(dipole_section, "REFERENCE_POINT", r_vals=t_control%dipole_ref_point)
1438 : ELSE
1439 6550 : NULLIFY (t_control%dipole_ref_point)
1440 6550 : IF (t_control%dipole_form == tddfpt_dipole_length .AND. t_control%dipole_reference == use_mom_ref_user) THEN
1441 0 : CPABORT("User-defined reference point should be given explicitly")
1442 : END IF
1443 : END IF
1444 :
1445 : !SOC subsection
1446 6550 : soc_section => section_vals_get_subs_vals(t_section, "SOC")
1447 6550 : CALL section_vals_get(soc_section, explicit=explicit)
1448 6550 : IF (explicit) THEN
1449 8 : t_control%do_soc = .TRUE.
1450 : END IF
1451 :
1452 : ! MGRID subsection
1453 6550 : mgrid_section => section_vals_get_subs_vals(t_section, "MGRID")
1454 6550 : CALL section_vals_get(mgrid_section, explicit=t_control%mgrid_is_explicit)
1455 :
1456 6550 : IF (t_control%mgrid_is_explicit) THEN
1457 2 : CALL section_vals_val_get(mgrid_section, "NGRIDS", i_val=t_control%mgrid_ngrids, explicit=explicit)
1458 2 : IF (.NOT. explicit) t_control%mgrid_ngrids = SIZE(qs_control%e_cutoff)
1459 :
1460 2 : CALL section_vals_val_get(mgrid_section, "CUTOFF", r_val=t_control%mgrid_cutoff, explicit=explicit)
1461 2 : IF (.NOT. explicit) t_control%mgrid_cutoff = qs_control%cutoff
1462 :
1463 : CALL section_vals_val_get(mgrid_section, "PROGRESSION_FACTOR", &
1464 2 : r_val=t_control%mgrid_progression_factor, explicit=explicit)
1465 2 : IF (explicit) THEN
1466 0 : IF (t_control%mgrid_progression_factor <= 1.0_dp) &
1467 : CALL cp_abort(__LOCATION__, &
1468 0 : "Progression factor should be greater then 1.0 to ensure multi-grid ordering")
1469 : ELSE
1470 2 : t_control%mgrid_progression_factor = qs_control%progression_factor
1471 : END IF
1472 :
1473 2 : CALL section_vals_val_get(mgrid_section, "COMMENSURATE", l_val=t_control%mgrid_commensurate_mgrids, explicit=explicit)
1474 2 : IF (.NOT. explicit) t_control%mgrid_commensurate_mgrids = qs_control%commensurate_mgrids
1475 2 : IF (t_control%mgrid_commensurate_mgrids) THEN
1476 0 : IF (explicit) THEN
1477 0 : t_control%mgrid_progression_factor = 4.0_dp
1478 : ELSE
1479 0 : t_control%mgrid_progression_factor = qs_control%progression_factor
1480 : END IF
1481 : END IF
1482 :
1483 2 : CALL section_vals_val_get(mgrid_section, "REL_CUTOFF", r_val=t_control%mgrid_relative_cutoff, explicit=explicit)
1484 2 : IF (.NOT. explicit) t_control%mgrid_relative_cutoff = qs_control%relative_cutoff
1485 :
1486 2 : CALL section_vals_val_get(mgrid_section, "MULTIGRID_SET", l_val=multigrid_set, explicit=explicit)
1487 2 : IF (.NOT. explicit) multigrid_set = .FALSE.
1488 2 : IF (multigrid_set) THEN
1489 0 : CALL section_vals_val_get(mgrid_section, "MULTIGRID_CUTOFF", r_vals=t_control%mgrid_e_cutoff)
1490 : ELSE
1491 2 : NULLIFY (t_control%mgrid_e_cutoff)
1492 : END IF
1493 :
1494 2 : CALL section_vals_val_get(mgrid_section, "REALSPACE", l_val=t_control%mgrid_realspace_mgrids, explicit=explicit)
1495 2 : IF (.NOT. explicit) t_control%mgrid_realspace_mgrids = qs_control%realspace_mgrids
1496 :
1497 : CALL section_vals_val_get(mgrid_section, "SKIP_LOAD_BALANCE_DISTRIBUTED", &
1498 2 : l_val=t_control%mgrid_skip_load_balance, explicit=explicit)
1499 2 : IF (.NOT. explicit) t_control%mgrid_skip_load_balance = qs_control%skip_load_balance_distributed
1500 :
1501 2 : IF (ASSOCIATED(t_control%mgrid_e_cutoff)) THEN
1502 0 : IF (SIZE(t_control%mgrid_e_cutoff) /= t_control%mgrid_ngrids) &
1503 0 : CPABORT("Inconsistent values for number of multi-grids")
1504 :
1505 : ! sort multi-grids in descending order according to their cutoff values
1506 0 : t_control%mgrid_e_cutoff = -t_control%mgrid_e_cutoff
1507 0 : ALLOCATE (inds(t_control%mgrid_ngrids))
1508 0 : CALL sort(t_control%mgrid_e_cutoff, t_control%mgrid_ngrids, inds)
1509 0 : DEALLOCATE (inds)
1510 0 : t_control%mgrid_e_cutoff = -t_control%mgrid_e_cutoff
1511 : END IF
1512 : END IF
1513 :
1514 : ! expand XC subsection (if given explicitly)
1515 6550 : xc_section => section_vals_get_subs_vals(t_section, "XC")
1516 6550 : xc_func => section_vals_get_subs_vals(xc_section, "XC_FUNCTIONAL")
1517 6550 : CALL section_vals_get(xc_func, explicit=explicit)
1518 6550 : IF (explicit) &
1519 192 : CALL xc_functionals_expand(xc_func, xc_section)
1520 :
1521 : ! sTDA subsection
1522 6550 : stda_section => section_vals_get_subs_vals(t_section, "STDA")
1523 6550 : IF (t_control%kernel == tddfpt_kernel_stda) THEN
1524 118 : t_control%stda_control%hfx_fraction = 0.0_dp
1525 118 : t_control%stda_control%do_exchange = .TRUE.
1526 118 : t_control%stda_control%eps_td_filter = 1.e-10_dp
1527 118 : t_control%stda_control%mn_alpha = -99.0_dp
1528 118 : t_control%stda_control%mn_beta = -99.0_dp
1529 : ! set default for Ewald method (on/off) dependent on periodicity
1530 212 : SELECT CASE (qs_control%periodicity)
1531 : CASE (0)
1532 94 : t_control%stda_control%do_ewald = .FALSE.
1533 : CASE (1)
1534 0 : t_control%stda_control%do_ewald = .TRUE.
1535 : CASE (2)
1536 0 : t_control%stda_control%do_ewald = .TRUE.
1537 : CASE (3)
1538 24 : t_control%stda_control%do_ewald = .TRUE.
1539 : CASE DEFAULT
1540 118 : CPABORT("Illegal value for periodiciy")
1541 : END SELECT
1542 118 : CALL section_vals_get(stda_section, explicit=explicit)
1543 118 : IF (explicit) THEN
1544 104 : CALL section_vals_val_get(stda_section, "HFX_FRACTION", r_val=hfx, explicit=expl)
1545 104 : IF (expl) t_control%stda_control%hfx_fraction = hfx
1546 104 : CALL section_vals_val_get(stda_section, "EPS_TD_FILTER", r_val=filter, explicit=expl)
1547 104 : IF (expl) t_control%stda_control%eps_td_filter = filter
1548 104 : CALL section_vals_val_get(stda_section, "DO_EWALD", l_val=do_ewald, explicit=expl)
1549 104 : IF (expl) t_control%stda_control%do_ewald = do_ewald
1550 104 : CALL section_vals_val_get(stda_section, "DO_EXCHANGE", l_val=do_exchange, explicit=expl)
1551 104 : IF (expl) t_control%stda_control%do_exchange = do_exchange
1552 104 : CALL section_vals_val_get(stda_section, "MATAGA_NISHIMOTO_CEXP", r_val=fval)
1553 104 : t_control%stda_control%mn_alpha = fval
1554 104 : CALL section_vals_val_get(stda_section, "MATAGA_NISHIMOTO_XEXP", r_val=fval)
1555 104 : t_control%stda_control%mn_beta = fval
1556 : END IF
1557 118 : CALL section_vals_val_get(stda_section, "COULOMB_SR_CUT", r_val=fval)
1558 118 : t_control%stda_control%coulomb_sr_cut = fval
1559 118 : CALL section_vals_val_get(stda_section, "COULOMB_SR_EPS", r_val=fval)
1560 118 : t_control%stda_control%coulomb_sr_eps = fval
1561 : END IF
1562 :
1563 6550 : CALL timestop(handle)
1564 6550 : END SUBROUTINE read_tddfpt2_control
1565 :
1566 : ! **************************************************************************************************
1567 : !> \brief Write the DFT control parameters to the output unit.
1568 : !> \param dft_control ...
1569 : !> \param dft_section ...
1570 : ! **************************************************************************************************
1571 11628 : SUBROUTINE write_dft_control(dft_control, dft_section)
1572 : TYPE(dft_control_type), POINTER :: dft_control
1573 : TYPE(section_vals_type), POINTER :: dft_section
1574 :
1575 : CHARACTER(len=*), PARAMETER :: routineN = 'write_dft_control'
1576 :
1577 : CHARACTER(LEN=20) :: tmpStr
1578 : INTEGER :: handle, output_unit
1579 : REAL(kind=dp) :: density_cut, density_smooth_cut_range, &
1580 : gradient_cut, tau_cut
1581 : TYPE(cp_logger_type), POINTER :: logger
1582 : TYPE(enumeration_type), POINTER :: enum
1583 : TYPE(keyword_type), POINTER :: keyword
1584 : TYPE(section_type), POINTER :: section
1585 : TYPE(section_vals_type), POINTER :: xc_section
1586 :
1587 7018 : IF (dft_control%qs_control%semi_empirical) RETURN
1588 5550 : IF (dft_control%qs_control%dftb) RETURN
1589 5328 : IF (dft_control%qs_control%xtb) THEN
1590 248 : CALL write_xtb_control(dft_control%qs_control%xtb_control, dft_section)
1591 248 : RETURN
1592 : END IF
1593 5080 : CALL timeset(routineN, handle)
1594 :
1595 5080 : NULLIFY (logger)
1596 5080 : logger => cp_get_default_logger()
1597 :
1598 : output_unit = cp_print_key_unit_nr(logger, dft_section, &
1599 5080 : "PRINT%DFT_CONTROL_PARAMETERS", extension=".Log")
1600 :
1601 5080 : IF (output_unit > 0) THEN
1602 :
1603 1300 : xc_section => section_vals_get_subs_vals(dft_section, "XC")
1604 :
1605 1300 : IF (dft_control%uks) THEN
1606 : WRITE (UNIT=output_unit, FMT="(/,T2,A,T78,A)") &
1607 406 : "DFT| Spin unrestricted (spin-polarized) Kohn-Sham calculation", "UKS"
1608 894 : ELSE IF (dft_control%roks) THEN
1609 : WRITE (UNIT=output_unit, FMT="(/,T2,A,T77,A)") &
1610 17 : "DFT| Spin restricted open Kohn-Sham calculation", "ROKS"
1611 : ELSE
1612 : WRITE (UNIT=output_unit, FMT="(/,T2,A,T78,A)") &
1613 877 : "DFT| Spin restricted Kohn-Sham (RKS) calculation", "RKS"
1614 : END IF
1615 :
1616 : WRITE (UNIT=output_unit, FMT="(T2,A,T76,I5)") &
1617 1300 : "DFT| Multiplicity", dft_control%multiplicity
1618 : WRITE (UNIT=output_unit, FMT="(T2,A,T76,I5)") &
1619 1300 : "DFT| Number of spin states", dft_control%nspins
1620 :
1621 : WRITE (UNIT=output_unit, FMT="(T2,A,T76,I5)") &
1622 1300 : "DFT| Charge", dft_control%charge
1623 :
1624 1300 : IF (dft_control%sic_method_id .NE. sic_none) CALL cite_reference(VandeVondele2005b)
1625 2586 : SELECT CASE (dft_control%sic_method_id)
1626 : CASE (sic_none)
1627 1286 : tmpstr = "NO"
1628 : CASE (sic_mauri_spz)
1629 6 : tmpstr = "SPZ/MAURI SIC"
1630 : CASE (sic_mauri_us)
1631 3 : tmpstr = "US/MAURI SIC"
1632 : CASE (sic_ad)
1633 3 : tmpstr = "AD SIC"
1634 : CASE (sic_eo)
1635 2 : tmpstr = "Explicit Orbital SIC"
1636 : CASE DEFAULT
1637 : ! fix throughout the cp2k for this option
1638 1300 : CPABORT("SIC option unknown")
1639 : END SELECT
1640 :
1641 : WRITE (UNIT=output_unit, FMT="(T2,A,T61,A20)") &
1642 1300 : "DFT| Self-interaction correction (SIC)", ADJUSTR(TRIM(tmpstr))
1643 :
1644 1300 : IF (dft_control%sic_method_id /= sic_none) THEN
1645 : WRITE (UNIT=output_unit, FMT="(T2,A,T66,ES15.6)") &
1646 14 : "DFT| SIC scaling parameter a", dft_control%sic_scaling_a, &
1647 28 : "DFT| SIC scaling parameter b", dft_control%sic_scaling_b
1648 : END IF
1649 :
1650 1300 : IF (dft_control%sic_method_id == sic_eo) THEN
1651 2 : IF (dft_control%sic_list_id == sic_list_all) THEN
1652 : WRITE (UNIT=output_unit, FMT="(T2,A,T66,A)") &
1653 1 : "DFT| SIC orbitals", "ALL"
1654 : END IF
1655 2 : IF (dft_control%sic_list_id == sic_list_unpaired) THEN
1656 : WRITE (UNIT=output_unit, FMT="(T2,A,T66,A)") &
1657 1 : "DFT| SIC orbitals", "UNPAIRED"
1658 : END IF
1659 : END IF
1660 :
1661 1300 : CALL section_vals_val_get(xc_section, "density_cutoff", r_val=density_cut)
1662 1300 : CALL section_vals_val_get(xc_section, "gradient_cutoff", r_val=gradient_cut)
1663 1300 : CALL section_vals_val_get(xc_section, "tau_cutoff", r_val=tau_cut)
1664 1300 : CALL section_vals_val_get(xc_section, "density_smooth_cutoff_range", r_val=density_smooth_cut_range)
1665 :
1666 : WRITE (UNIT=output_unit, FMT="(T2,A,T66,ES15.6)") &
1667 1300 : "DFT| Cutoffs: density ", density_cut, &
1668 1300 : "DFT| gradient", gradient_cut, &
1669 1300 : "DFT| tau ", tau_cut, &
1670 2600 : "DFT| cutoff_smoothing_range", density_smooth_cut_range
1671 : CALL section_vals_val_get(xc_section, "XC_GRID%XC_SMOOTH_RHO", &
1672 1300 : c_val=tmpStr)
1673 : WRITE (output_unit, '( A, T61, A )') &
1674 1300 : " DFT| XC density smoothing ", ADJUSTR(tmpStr)
1675 : CALL section_vals_val_get(xc_section, "XC_GRID%XC_DERIV", &
1676 1300 : c_val=tmpStr)
1677 : WRITE (output_unit, '( A, T61, A )') &
1678 1300 : " DFT| XC derivatives ", ADJUSTR(tmpStr)
1679 1300 : IF (dft_control%dft_plus_u) THEN
1680 16 : NULLIFY (enum, keyword, section)
1681 16 : CALL create_dft_section(section)
1682 16 : keyword => section_get_keyword(section, "PLUS_U_METHOD")
1683 16 : CALL keyword_get(keyword, enum=enum)
1684 : WRITE (UNIT=output_unit, FMT="(/,T2,A,T41,A40)") &
1685 16 : "DFT+U| Method", ADJUSTR(TRIM(enum_i2c(enum, dft_control%plus_u_method_id)))
1686 : WRITE (UNIT=output_unit, FMT="(T2,A)") &
1687 16 : "DFT+U| Check atomic kind information for details"
1688 16 : CALL section_release(section)
1689 : END IF
1690 :
1691 1300 : WRITE (UNIT=output_unit, FMT="(A)") ""
1692 1300 : CALL xc_write(output_unit, xc_section, dft_control%lsd)
1693 :
1694 1300 : IF (dft_control%apply_period_efield) THEN
1695 4 : WRITE (UNIT=output_unit, FMT="(A)") ""
1696 4 : IF (dft_control%period_efield%displacement_field) THEN
1697 : WRITE (UNIT=output_unit, FMT="(T2,A)") &
1698 0 : "PERIODIC_EFIELD| Use displacement field formulation"
1699 : WRITE (UNIT=output_unit, FMT="(T2,A,T66,1X,ES14.6)") &
1700 0 : "PERIODIC_EFIELD| Displacement field filter: x", &
1701 0 : dft_control%period_efield%d_filter(1), &
1702 0 : "PERIODIC_EFIELD| y", &
1703 0 : dft_control%period_efield%d_filter(2), &
1704 0 : "PERIODIC_EFIELD| z", &
1705 0 : dft_control%period_efield%d_filter(3)
1706 : END IF
1707 : WRITE (UNIT=output_unit, FMT="(T2,A,T66,1X,ES14.6)") &
1708 4 : "PERIODIC_EFIELD| Polarisation vector: x", &
1709 4 : dft_control%period_efield%polarisation(1), &
1710 4 : "PERIODIC_EFIELD| y", &
1711 4 : dft_control%period_efield%polarisation(2), &
1712 4 : "PERIODIC_EFIELD| z", &
1713 4 : dft_control%period_efield%polarisation(3), &
1714 4 : "PERIODIC_EFIELD| Intensity [a.u.]:", &
1715 8 : dft_control%period_efield%strength
1716 16 : IF (SQRT(DOT_PRODUCT(dft_control%period_efield%polarisation, &
1717 : dft_control%period_efield%polarisation)) < EPSILON(0.0_dp)) THEN
1718 0 : CPABORT("Invalid (too small) polarisation vector specified for PERIODIC_EFIELD")
1719 : END IF
1720 : END IF
1721 :
1722 1300 : IF (dft_control%do_sccs) THEN
1723 : WRITE (UNIT=output_unit, FMT="(/,T2,A)") &
1724 5 : "SCCS| Self-consistent continuum solvation model"
1725 : WRITE (UNIT=output_unit, FMT="(T2,A,T61,ES20.6)") &
1726 5 : "SCCS| Relative permittivity of the solvent (medium)", &
1727 5 : dft_control%sccs_control%epsilon_solvent, &
1728 5 : "SCCS| Absolute permittivity [a.u.]", &
1729 10 : dft_control%sccs_control%epsilon_solvent/fourpi
1730 9 : SELECT CASE (dft_control%sccs_control%method_id)
1731 : CASE (sccs_andreussi)
1732 : WRITE (UNIT=output_unit, FMT="(T2,A,/,(T2,A,T61,ES20.6))") &
1733 4 : "SCCS| Dielectric function proposed by Andreussi et al.", &
1734 4 : "SCCS| rho_max", dft_control%sccs_control%rho_max, &
1735 8 : "SCCS| rho_min", dft_control%sccs_control%rho_min
1736 : CASE (sccs_fattebert_gygi)
1737 : WRITE (UNIT=output_unit, FMT="(T2,A,/,(T2,A,T61,ES20.6))") &
1738 1 : "SCCS| Dielectric function proposed by Fattebert and Gygi", &
1739 1 : "SCCS| beta", dft_control%sccs_control%beta, &
1740 2 : "SCCS| rho_zero", dft_control%sccs_control%rho_zero
1741 : CASE DEFAULT
1742 5 : CPABORT("Invalid SCCS model specified. Please, check your input!")
1743 : END SELECT
1744 6 : SELECT CASE (dft_control%sccs_control%derivative_method)
1745 : CASE (sccs_derivative_fft)
1746 : WRITE (UNIT=output_unit, FMT="(T2,A,T46,A35)") &
1747 1 : "SCCS| Numerical derivative calculation", &
1748 2 : ADJUSTR("FFT")
1749 : CASE (sccs_derivative_cd3)
1750 : WRITE (UNIT=output_unit, FMT="(T2,A,T46,A35)") &
1751 0 : "SCCS| Numerical derivative calculation", &
1752 0 : ADJUSTR("3-point stencil central differences")
1753 : CASE (sccs_derivative_cd5)
1754 : WRITE (UNIT=output_unit, FMT="(T2,A,T46,A35)") &
1755 4 : "SCCS| Numerical derivative calculation", &
1756 8 : ADJUSTR("5-point stencil central differences")
1757 : CASE (sccs_derivative_cd7)
1758 : WRITE (UNIT=output_unit, FMT="(T2,A,T46,A35)") &
1759 0 : "SCCS| Numerical derivative calculation", &
1760 0 : ADJUSTR("7-point stencil central differences")
1761 : CASE DEFAULT
1762 : CALL cp_abort(__LOCATION__, &
1763 : "Invalid derivative method specified for SCCS model. "// &
1764 5 : "Please, check your input!")
1765 : END SELECT
1766 : WRITE (UNIT=output_unit, FMT="(T2,A,T61,ES20.6)") &
1767 5 : "SCCS| Repulsion parameter alpha [mN/m] = [dyn/cm]", &
1768 10 : cp_unit_from_cp2k(dft_control%sccs_control%alpha_solvent, "mN/m")
1769 : WRITE (UNIT=output_unit, FMT="(T2,A,T61,ES20.6)") &
1770 5 : "SCCS| Dispersion parameter beta [GPa]", &
1771 10 : cp_unit_from_cp2k(dft_control%sccs_control%beta_solvent, "GPa")
1772 : WRITE (UNIT=output_unit, FMT="(T2,A,T61,ES20.6)") &
1773 5 : "SCCS| Surface tension gamma [mN/m] = [dyn/cm]", &
1774 10 : cp_unit_from_cp2k(dft_control%sccs_control%gamma_solvent, "mN/m")
1775 : WRITE (UNIT=output_unit, FMT="(T2,A,T61,ES20.6)") &
1776 5 : "SCCS| Mixing parameter applied during the iteration cycle", &
1777 10 : dft_control%sccs_control%mixing
1778 : WRITE (UNIT=output_unit, FMT="(T2,A,T61,ES20.6)") &
1779 5 : "SCCS| Tolerance for the convergence of the SCCS iteration cycle", &
1780 10 : dft_control%sccs_control%eps_sccs
1781 : WRITE (UNIT=output_unit, FMT="(T2,A,T61,I20)") &
1782 5 : "SCCS| Maximum number of iteration steps", &
1783 10 : dft_control%sccs_control%max_iter
1784 : WRITE (UNIT=output_unit, FMT="(T2,A,T61,ES20.6)") &
1785 5 : "SCCS| SCF convergence threshold for starting the SCCS iteration", &
1786 10 : dft_control%sccs_control%eps_scf
1787 : WRITE (UNIT=output_unit, FMT="(T2,A,T61,ES20.6)") &
1788 5 : "SCCS| Numerical increment for the cavity surface calculation", &
1789 10 : dft_control%sccs_control%delta_rho
1790 : END IF
1791 :
1792 1300 : WRITE (UNIT=output_unit, FMT="(A)") ""
1793 :
1794 : END IF
1795 :
1796 : CALL cp_print_key_finished_output(output_unit, logger, dft_section, &
1797 5080 : "PRINT%DFT_CONTROL_PARAMETERS")
1798 :
1799 5080 : CALL timestop(handle)
1800 :
1801 : END SUBROUTINE write_dft_control
1802 :
1803 : ! **************************************************************************************************
1804 : !> \brief Write the ADMM control parameters to the output unit.
1805 : !> \param admm_control ...
1806 : !> \param dft_section ...
1807 : ! **************************************************************************************************
1808 442 : SUBROUTINE write_admm_control(admm_control, dft_section)
1809 : TYPE(admm_control_type), POINTER :: admm_control
1810 : TYPE(section_vals_type), POINTER :: dft_section
1811 :
1812 : INTEGER :: iounit
1813 : TYPE(cp_logger_type), POINTER :: logger
1814 :
1815 442 : NULLIFY (logger)
1816 442 : logger => cp_get_default_logger()
1817 :
1818 : iounit = cp_print_key_unit_nr(logger, dft_section, &
1819 442 : "PRINT%DFT_CONTROL_PARAMETERS", extension=".Log")
1820 :
1821 442 : IF (iounit > 0) THEN
1822 :
1823 233 : SELECT CASE (admm_control%admm_type)
1824 : CASE (no_admm_type)
1825 114 : WRITE (UNIT=iounit, FMT="(/,T2,A,T77,A)") "ADMM| Specific ADMM type specified", "NONE"
1826 : CASE (admm1_type)
1827 1 : WRITE (UNIT=iounit, FMT="(/,T2,A,T76,A)") "ADMM| Specific ADMM type specified", "ADMM1"
1828 : CASE (admm2_type)
1829 1 : WRITE (UNIT=iounit, FMT="(/,T2,A,T76,A)") "ADMM| Specific ADMM type specified", "ADMM2"
1830 : CASE (admms_type)
1831 1 : WRITE (UNIT=iounit, FMT="(/,T2,A,T76,A)") "ADMM| Specific ADMM type specified", "ADMMS"
1832 : CASE (admmp_type)
1833 1 : WRITE (UNIT=iounit, FMT="(/,T2,A,T76,A)") "ADMM| Specific ADMM type specified", "ADMMP"
1834 : CASE (admmq_type)
1835 1 : WRITE (UNIT=iounit, FMT="(/,T2,A,T76,A)") "ADMM| Specific ADMM type specified", "ADMMQ"
1836 : CASE DEFAULT
1837 119 : CPABORT("admm_type")
1838 : END SELECT
1839 :
1840 189 : SELECT CASE (admm_control%purification_method)
1841 : CASE (do_admm_purify_none)
1842 70 : WRITE (UNIT=iounit, FMT="(T2,A,T77,A)") "ADMM| Density matrix purification method", "NONE"
1843 : CASE (do_admm_purify_cauchy)
1844 9 : WRITE (UNIT=iounit, FMT="(T2,A,T75,A)") "ADMM| Density matrix purification method", "Cauchy"
1845 : CASE (do_admm_purify_cauchy_subspace)
1846 5 : WRITE (UNIT=iounit, FMT="(T2,A,T66,A)") "ADMM| Density matrix purification method", "Cauchy subspace"
1847 : CASE (do_admm_purify_mo_diag)
1848 25 : WRITE (UNIT=iounit, FMT="(T2,A,T63,A)") "ADMM| Density matrix purification method", "MO diagonalization"
1849 : CASE (do_admm_purify_mo_no_diag)
1850 3 : WRITE (UNIT=iounit, FMT="(T2,A,T71,A)") "ADMM| Density matrix purification method", "MO no diag"
1851 : CASE (do_admm_purify_mcweeny)
1852 1 : WRITE (UNIT=iounit, FMT="(T2,A,T74,A)") "ADMM| Density matrix purification method", "McWeeny"
1853 : CASE (do_admm_purify_none_dm)
1854 6 : WRITE (UNIT=iounit, FMT="(T2,A,T73,A)") "ADMM| Density matrix purification method", "NONE(DM)"
1855 : CASE DEFAULT
1856 119 : CPABORT("admm_purification_method")
1857 : END SELECT
1858 :
1859 214 : SELECT CASE (admm_control%method)
1860 : CASE (do_admm_basis_projection)
1861 95 : WRITE (UNIT=iounit, FMT="(T2,A)") "ADMM| Orbital projection on ADMM basis"
1862 : CASE (do_admm_blocking_purify_full)
1863 3 : WRITE (UNIT=iounit, FMT="(T2,A)") "ADMM| Blocked Fock matrix projection with full purification"
1864 : CASE (do_admm_blocked_projection)
1865 6 : WRITE (UNIT=iounit, FMT="(T2,A)") "ADMM| Blocked Fock matrix projection"
1866 : CASE (do_admm_charge_constrained_projection)
1867 15 : WRITE (UNIT=iounit, FMT="(T2,A)") "ADMM| Orbital projection with charge constrain"
1868 : CASE DEFAULT
1869 119 : CPABORT("admm method")
1870 : END SELECT
1871 :
1872 136 : SELECT CASE (admm_control%scaling_model)
1873 : CASE (do_admm_exch_scaling_none)
1874 : CASE (do_admm_exch_scaling_merlot)
1875 17 : WRITE (UNIT=iounit, FMT="(T2,A)") "ADMM| Use Merlot (2014) scaling model"
1876 : CASE DEFAULT
1877 119 : CPABORT("admm scaling_model")
1878 : END SELECT
1879 :
1880 119 : WRITE (UNIT=iounit, FMT="(T2,A,T61,G20.10)") "ADMM| eps_filter", admm_control%eps_filter
1881 :
1882 127 : SELECT CASE (admm_control%aux_exch_func)
1883 : CASE (do_admm_aux_exch_func_none)
1884 8 : WRITE (UNIT=iounit, FMT="(T2,A)") "ADMM| No exchange functional correction term used"
1885 : CASE (do_admm_aux_exch_func_default, do_admm_aux_exch_func_default_libxc)
1886 85 : WRITE (UNIT=iounit, FMT="(T2,A,T74,A)") "ADMM| Exchange functional in correction term", "(W)PBEX"
1887 : CASE (do_admm_aux_exch_func_pbex, do_admm_aux_exch_func_pbex_libxc)
1888 17 : WRITE (UNIT=iounit, FMT="(T2,A,T77,A)") "ADMM| Exchange functional in correction term", "PBEX"
1889 : CASE (do_admm_aux_exch_func_opt, do_admm_aux_exch_func_opt_libxc)
1890 8 : WRITE (UNIT=iounit, FMT="(T2,A,T77,A)") "ADMM| Exchange functional in correction term", "OPTX"
1891 : CASE (do_admm_aux_exch_func_bee, do_admm_aux_exch_func_bee_libxc)
1892 1 : WRITE (UNIT=iounit, FMT="(T2,A,T74,A)") "ADMM| Exchange functional in correction term", "Becke88"
1893 : CASE (do_admm_aux_exch_func_sx_libxc)
1894 0 : WRITE (UNIT=iounit, FMT="(T2,A,T74,A)") "ADMM| Exchange functional in correction term", "SlaterX"
1895 : CASE DEFAULT
1896 119 : CPABORT("admm aux_exch_func")
1897 : END SELECT
1898 :
1899 119 : WRITE (UNIT=iounit, FMT="(A)") ""
1900 :
1901 : END IF
1902 :
1903 : CALL cp_print_key_finished_output(iounit, logger, dft_section, &
1904 442 : "PRINT%DFT_CONTROL_PARAMETERS")
1905 442 : END SUBROUTINE write_admm_control
1906 :
1907 : ! **************************************************************************************************
1908 : !> \brief Write the xTB control parameters to the output unit.
1909 : !> \param xtb_control ...
1910 : !> \param dft_section ...
1911 : ! **************************************************************************************************
1912 248 : SUBROUTINE write_xtb_control(xtb_control, dft_section)
1913 : TYPE(xtb_control_type), POINTER :: xtb_control
1914 : TYPE(section_vals_type), POINTER :: dft_section
1915 :
1916 : CHARACTER(len=*), PARAMETER :: routineN = 'write_xtb_control'
1917 :
1918 : INTEGER :: handle, output_unit
1919 : TYPE(cp_logger_type), POINTER :: logger
1920 :
1921 248 : CALL timeset(routineN, handle)
1922 248 : NULLIFY (logger)
1923 248 : logger => cp_get_default_logger()
1924 :
1925 : output_unit = cp_print_key_unit_nr(logger, dft_section, &
1926 248 : "PRINT%DFT_CONTROL_PARAMETERS", extension=".Log")
1927 :
1928 248 : IF (output_unit > 0) THEN
1929 :
1930 : WRITE (UNIT=output_unit, FMT="(/,T2,A,T31,A50)") &
1931 37 : "xTB| Parameter file", ADJUSTR(TRIM(xtb_control%parameter_file_name))
1932 : WRITE (UNIT=output_unit, FMT="(T2,A,T71,I10)") &
1933 37 : "xTB| Basis expansion STO-NG", xtb_control%sto_ng
1934 : WRITE (UNIT=output_unit, FMT="(T2,A,T71,I10)") &
1935 37 : "xTB| Basis expansion STO-NG for Hydrogen", xtb_control%h_sto_ng
1936 : WRITE (UNIT=output_unit, FMT="(T2,A,T71,L10)") &
1937 37 : "xTB| Halogen interaction potential", xtb_control%xb_interaction
1938 : WRITE (UNIT=output_unit, FMT="(T2,A,T71,F10.3)") &
1939 37 : "xTB| Halogen interaction potential cutoff radius", xtb_control%xb_radius
1940 : WRITE (UNIT=output_unit, FMT="(T2,A,T71,L10)") &
1941 37 : "xTB| Nonbonded interactions", xtb_control%do_nonbonded
1942 : WRITE (UNIT=output_unit, FMT="(T2,A,T31,A50)") &
1943 37 : "xTB| D3 Dispersion: Parameter file", ADJUSTR(TRIM(xtb_control%dispersion_parameter_file))
1944 : WRITE (UNIT=output_unit, FMT="(T2,A,T51,3F10.3)") &
1945 37 : "xTB| Huckel constants ks kp kd", xtb_control%ks, xtb_control%kp, xtb_control%kd
1946 : WRITE (UNIT=output_unit, FMT="(T2,A,T61,2F10.3)") &
1947 37 : "xTB| Huckel constants ksp k2sh", xtb_control%ksp, xtb_control%k2sh
1948 : WRITE (UNIT=output_unit, FMT="(T2,A,T71,F10.3)") &
1949 37 : "xTB| Mataga-Nishimoto exponent", xtb_control%kg
1950 : WRITE (UNIT=output_unit, FMT="(T2,A,T71,F10.3)") &
1951 37 : "xTB| Repulsion potential exponent", xtb_control%kf
1952 : WRITE (UNIT=output_unit, FMT="(T2,A,T51,3F10.3)") &
1953 37 : "xTB| Coordination number scaling kcn(s) kcn(p) kcn(d)", &
1954 74 : xtb_control%kcns, xtb_control%kcnp, xtb_control%kcnd
1955 : WRITE (UNIT=output_unit, FMT="(T2,A,T71,F10.3)") &
1956 37 : "xTB| Electronegativity scaling", xtb_control%ken
1957 : WRITE (UNIT=output_unit, FMT="(T2,A,T61,2F10.3)") &
1958 37 : "xTB| Halogen potential scaling kxr kx2", xtb_control%kxr, xtb_control%kx2
1959 37 : WRITE (UNIT=output_unit, FMT="(/)")
1960 :
1961 : END IF
1962 :
1963 : CALL cp_print_key_finished_output(output_unit, logger, dft_section, &
1964 248 : "PRINT%DFT_CONTROL_PARAMETERS")
1965 :
1966 248 : CALL timestop(handle)
1967 :
1968 248 : END SUBROUTINE write_xtb_control
1969 :
1970 : ! **************************************************************************************************
1971 : !> \brief Purpose: Write the QS control parameters to the output unit.
1972 : !> \param qs_control ...
1973 : !> \param dft_section ...
1974 : ! **************************************************************************************************
1975 11628 : SUBROUTINE write_qs_control(qs_control, dft_section)
1976 : TYPE(qs_control_type), INTENT(IN) :: qs_control
1977 : TYPE(section_vals_type), POINTER :: dft_section
1978 :
1979 : CHARACTER(len=*), PARAMETER :: routineN = 'write_qs_control'
1980 :
1981 : CHARACTER(len=20) :: method, quadrature
1982 : INTEGER :: handle, i, igrid_level, ngrid_level, &
1983 : output_unit
1984 : TYPE(cp_logger_type), POINTER :: logger
1985 : TYPE(ddapc_restraint_type), POINTER :: ddapc_restraint_control
1986 : TYPE(enumeration_type), POINTER :: enum
1987 : TYPE(keyword_type), POINTER :: keyword
1988 : TYPE(section_type), POINTER :: qs_section
1989 : TYPE(section_vals_type), POINTER :: print_section_vals, qs_section_vals
1990 :
1991 7018 : IF (qs_control%semi_empirical) RETURN
1992 5550 : IF (qs_control%dftb) RETURN
1993 5328 : IF (qs_control%xtb) RETURN
1994 5080 : CALL timeset(routineN, handle)
1995 5080 : NULLIFY (logger, print_section_vals, qs_section, qs_section_vals)
1996 5080 : logger => cp_get_default_logger()
1997 5080 : print_section_vals => section_vals_get_subs_vals(dft_section, "PRINT")
1998 5080 : qs_section_vals => section_vals_get_subs_vals(dft_section, "QS")
1999 5080 : CALL section_vals_get(qs_section_vals, section=qs_section)
2000 :
2001 5080 : NULLIFY (enum, keyword)
2002 5080 : keyword => section_get_keyword(qs_section, "METHOD")
2003 5080 : CALL keyword_get(keyword, enum=enum)
2004 5080 : method = enum_i2c(enum, qs_control%method_id)
2005 :
2006 5080 : NULLIFY (enum, keyword)
2007 5080 : keyword => section_get_keyword(qs_section, "QUADRATURE")
2008 5080 : CALL keyword_get(keyword, enum=enum)
2009 5080 : quadrature = enum_i2c(enum, qs_control%gapw_control%quadrature)
2010 :
2011 : output_unit = cp_print_key_unit_nr(logger, print_section_vals, &
2012 5080 : "DFT_CONTROL_PARAMETERS", extension=".Log")
2013 5080 : IF (output_unit > 0) THEN
2014 1300 : ngrid_level = SIZE(qs_control%e_cutoff)
2015 : WRITE (UNIT=output_unit, FMT="(/,T2,A,T61,A20)") &
2016 1300 : "QS| Method:", ADJUSTR(method)
2017 1300 : IF (qs_control%pw_grid_opt%spherical) THEN
2018 : WRITE (UNIT=output_unit, FMT="(T2,A,T61,A)") &
2019 0 : "QS| Density plane wave grid type", " SPHERICAL HALFSPACE"
2020 1300 : ELSE IF (qs_control%pw_grid_opt%fullspace) THEN
2021 : WRITE (UNIT=output_unit, FMT="(T2,A,T57,A)") &
2022 1300 : "QS| Density plane wave grid type", " NON-SPHERICAL FULLSPACE"
2023 : ELSE
2024 : WRITE (UNIT=output_unit, FMT="(T2,A,T57,A)") &
2025 0 : "QS| Density plane wave grid type", " NON-SPHERICAL HALFSPACE"
2026 : END IF
2027 : WRITE (UNIT=output_unit, FMT="(T2,A,T71,I10)") &
2028 1300 : "QS| Number of grid levels:", SIZE(qs_control%e_cutoff)
2029 1300 : IF (ngrid_level == 1) THEN
2030 : WRITE (UNIT=output_unit, FMT="(T2,A,T71,F10.1)") &
2031 61 : "QS| Density cutoff [a.u.]:", qs_control%e_cutoff(1)
2032 : ELSE
2033 : WRITE (UNIT=output_unit, FMT="(T2,A,T71,F10.1)") &
2034 1239 : "QS| Density cutoff [a.u.]:", qs_control%cutoff
2035 1239 : IF (qs_control%commensurate_mgrids) &
2036 131 : WRITE (UNIT=output_unit, FMT="(T2,A)") "QS| Using commensurate multigrids"
2037 : WRITE (UNIT=output_unit, FMT="(T2,A,T71,F10.1)") &
2038 1239 : "QS| Multi grid cutoff [a.u.]: 1) grid level", qs_control%e_cutoff(1)
2039 : WRITE (UNIT=output_unit, FMT="(T2,A,I3,A,T71,F10.1)") &
2040 3882 : ("QS| ", igrid_level, ") grid level", &
2041 5121 : qs_control%e_cutoff(igrid_level), &
2042 6360 : igrid_level=2, SIZE(qs_control%e_cutoff))
2043 : END IF
2044 1300 : IF (qs_control%pao) THEN
2045 0 : WRITE (UNIT=output_unit, FMT="(T2,A)") "QS| PAO active"
2046 : END IF
2047 : WRITE (UNIT=output_unit, FMT="(T2,A,T71,F10.1)") &
2048 1300 : "QS| Grid level progression factor:", qs_control%progression_factor
2049 : WRITE (UNIT=output_unit, FMT="(T2,A,T71,F10.1)") &
2050 1300 : "QS| Relative density cutoff [a.u.]:", qs_control%relative_cutoff
2051 : WRITE (UNIT=output_unit, FMT="(T2,A,T73,ES8.1)") &
2052 1300 : "QS| Interaction thresholds: eps_pgf_orb:", &
2053 1300 : qs_control%eps_pgf_orb, &
2054 1300 : "QS| eps_filter_matrix:", &
2055 1300 : qs_control%eps_filter_matrix, &
2056 1300 : "QS| eps_core_charge:", &
2057 1300 : qs_control%eps_core_charge, &
2058 1300 : "QS| eps_rho_gspace:", &
2059 1300 : qs_control%eps_rho_gspace, &
2060 1300 : "QS| eps_rho_rspace:", &
2061 1300 : qs_control%eps_rho_rspace, &
2062 1300 : "QS| eps_gvg_rspace:", &
2063 1300 : qs_control%eps_gvg_rspace, &
2064 1300 : "QS| eps_ppl:", &
2065 1300 : qs_control%eps_ppl, &
2066 1300 : "QS| eps_ppnl:", &
2067 2600 : qs_control%eps_ppnl
2068 1300 : IF (qs_control%gapw) THEN
2069 380 : SELECT CASE (qs_control%gapw_control%basis_1c)
2070 : CASE (gapw_1c_orb)
2071 : WRITE (UNIT=output_unit, FMT="(T2,A)") &
2072 187 : "QS| GAPW| One center basis from orbital basis primitives"
2073 : CASE (gapw_1c_small)
2074 : WRITE (UNIT=output_unit, FMT="(T2,A)") &
2075 2 : "QS| GAPW| One center basis extended with primitives (small:s)"
2076 : CASE (gapw_1c_medium)
2077 : WRITE (UNIT=output_unit, FMT="(T2,A)") &
2078 1 : "QS| GAPW| One center basis extended with primitives (medium:sp)"
2079 : CASE (gapw_1c_large)
2080 : WRITE (UNIT=output_unit, FMT="(T2,A)") &
2081 2 : "QS| GAPW| One center basis extended with primitives (large:spd)"
2082 : CASE (gapw_1c_very_large)
2083 : WRITE (UNIT=output_unit, FMT="(T2,A)") &
2084 1 : "QS| GAPW| One center basis extended with primitives (very large:spdf)"
2085 : CASE DEFAULT
2086 193 : CPABORT("basis_1c incorrect")
2087 : END SELECT
2088 : WRITE (UNIT=output_unit, FMT="(T2,A,T73,ES8.1)") &
2089 193 : "QS| GAPW| eps_fit:", &
2090 193 : qs_control%gapw_control%eps_fit, &
2091 193 : "QS| GAPW| eps_iso:", &
2092 193 : qs_control%gapw_control%eps_iso, &
2093 193 : "QS| GAPW| eps_svd:", &
2094 193 : qs_control%gapw_control%eps_svd, &
2095 193 : "QS| GAPW| eps_cpc:", &
2096 386 : qs_control%gapw_control%eps_cpc
2097 : WRITE (UNIT=output_unit, FMT="(T2,A,T61,A20)") &
2098 193 : "QS| GAPW| atom-r-grid: quadrature:", &
2099 386 : ADJUSTR(quadrature)
2100 : WRITE (UNIT=output_unit, FMT="(T2,A,T71,I10)") &
2101 193 : "QS| GAPW| atom-s-grid: max l :", &
2102 193 : qs_control%gapw_control%lmax_sphere, &
2103 193 : "QS| GAPW| max_l_rho0 :", &
2104 386 : qs_control%gapw_control%lmax_rho0
2105 193 : IF (qs_control%gapw_control%non_paw_atoms) THEN
2106 : WRITE (UNIT=output_unit, FMT="(T2,A)") &
2107 29 : "QS| GAPW| At least one kind is NOT PAW, i.e. it has only soft AO "
2108 : END IF
2109 193 : IF (qs_control%gapw_control%nopaw_as_gpw) THEN
2110 : WRITE (UNIT=output_unit, FMT="(T2,A)") &
2111 29 : "QS| GAPW| The NOT PAW atoms are treated fully GPW"
2112 : END IF
2113 : END IF
2114 1300 : IF (qs_control%gapw_xc) THEN
2115 50 : SELECT CASE (qs_control%gapw_control%basis_1c)
2116 : CASE (gapw_1c_orb)
2117 : WRITE (UNIT=output_unit, FMT="(T2,A)") &
2118 25 : "QS| GAPW_XC| One center basis from orbital basis primitives"
2119 : CASE (gapw_1c_small)
2120 : WRITE (UNIT=output_unit, FMT="(T2,A)") &
2121 0 : "QS| GAPW_XC| One center basis extended with primitives (small:s)"
2122 : CASE (gapw_1c_medium)
2123 : WRITE (UNIT=output_unit, FMT="(T2,A)") &
2124 0 : "QS| GAPW_XC| One center basis extended with primitives (medium:sp)"
2125 : CASE (gapw_1c_large)
2126 : WRITE (UNIT=output_unit, FMT="(T2,A)") &
2127 0 : "QS| GAPW_XC| One center basis extended with primitives (large:spd)"
2128 : CASE (gapw_1c_very_large)
2129 : WRITE (UNIT=output_unit, FMT="(T2,A)") &
2130 0 : "QS| GAPW_XC| One center basis extended with primitives (very large:spdf)"
2131 : CASE DEFAULT
2132 25 : CPABORT("basis_1c incorrect")
2133 : END SELECT
2134 : WRITE (UNIT=output_unit, FMT="(T2,A,T73,ES8.1)") &
2135 25 : "QS| GAPW_XC| eps_fit:", &
2136 25 : qs_control%gapw_control%eps_fit, &
2137 25 : "QS| GAPW_XC| eps_iso:", &
2138 25 : qs_control%gapw_control%eps_iso, &
2139 25 : "QS| GAPW_XC| eps_svd:", &
2140 50 : qs_control%gapw_control%eps_svd
2141 : WRITE (UNIT=output_unit, FMT="(T2,A,T55,A30)") &
2142 25 : "QS| GAPW_XC|atom-r-grid: quadrature:", &
2143 50 : enum_i2c(enum, qs_control%gapw_control%quadrature)
2144 : WRITE (UNIT=output_unit, FMT="(T2,A,T71,I10)") &
2145 25 : "QS| GAPW_XC| atom-s-grid: max l :", &
2146 50 : qs_control%gapw_control%lmax_sphere
2147 : END IF
2148 1300 : IF (qs_control%mulliken_restraint) THEN
2149 : WRITE (UNIT=output_unit, FMT="(T2,A,T73,ES8.1)") &
2150 1 : "QS| Mulliken restraint target", qs_control%mulliken_restraint_control%target
2151 : WRITE (UNIT=output_unit, FMT="(T2,A,T73,ES8.1)") &
2152 1 : "QS| Mulliken restraint strength", qs_control%mulliken_restraint_control%strength
2153 : WRITE (UNIT=output_unit, FMT="(T2,A,T73,I8)") &
2154 1 : "QS| Mulliken restraint atoms: ", qs_control%mulliken_restraint_control%natoms
2155 2 : WRITE (UNIT=output_unit, FMT="(5I8)") qs_control%mulliken_restraint_control%atoms
2156 : END IF
2157 1300 : IF (qs_control%ddapc_restraint) THEN
2158 14 : DO i = 1, SIZE(qs_control%ddapc_restraint_control)
2159 8 : ddapc_restraint_control => qs_control%ddapc_restraint_control(i)
2160 8 : IF (SIZE(qs_control%ddapc_restraint_control) .GT. 1) &
2161 : WRITE (UNIT=output_unit, FMT="(T2,A,T3,I8)") &
2162 3 : "QS| parameters for DDAPC restraint number", i
2163 : WRITE (UNIT=output_unit, FMT="(T2,A,T73,ES8.1)") &
2164 8 : "QS| ddapc restraint target", ddapc_restraint_control%target
2165 : WRITE (UNIT=output_unit, FMT="(T2,A,T73,ES8.1)") &
2166 8 : "QS| ddapc restraint strength", ddapc_restraint_control%strength
2167 : WRITE (UNIT=output_unit, FMT="(T2,A,T73,I8)") &
2168 8 : "QS| ddapc restraint atoms: ", ddapc_restraint_control%natoms
2169 17 : WRITE (UNIT=output_unit, FMT="(5I8)") ddapc_restraint_control%atoms
2170 8 : WRITE (UNIT=output_unit, FMT="(T2,A)") "Coefficients:"
2171 17 : WRITE (UNIT=output_unit, FMT="(5F6.2)") ddapc_restraint_control%coeff
2172 6 : SELECT CASE (ddapc_restraint_control%functional_form)
2173 : CASE (do_ddapc_restraint)
2174 : WRITE (UNIT=output_unit, FMT="(T2,A,T61,A20)") &
2175 3 : "QS| ddapc restraint functional form :", "RESTRAINT"
2176 : CASE (do_ddapc_constraint)
2177 : WRITE (UNIT=output_unit, FMT="(T2,A,T61,A20)") &
2178 5 : "QS| ddapc restraint functional form :", "CONSTRAINT"
2179 : CASE DEFAULT
2180 8 : CPABORT("Unknown ddapc restraint")
2181 : END SELECT
2182 : END DO
2183 : END IF
2184 1300 : IF (qs_control%s2_restraint) THEN
2185 : WRITE (UNIT=output_unit, FMT="(T2,A,T73,ES8.1)") &
2186 0 : "QS| s2 restraint target", qs_control%s2_restraint_control%target
2187 : WRITE (UNIT=output_unit, FMT="(T2,A,T73,ES8.1)") &
2188 0 : "QS| s2 restraint strength", qs_control%s2_restraint_control%strength
2189 0 : SELECT CASE (qs_control%s2_restraint_control%functional_form)
2190 : CASE (do_s2_restraint)
2191 : WRITE (UNIT=output_unit, FMT="(T2,A,T61,A20)") &
2192 0 : "QS| s2 restraint functional form :", "RESTRAINT"
2193 0 : CPABORT("Not yet implemented")
2194 : CASE (do_s2_constraint)
2195 : WRITE (UNIT=output_unit, FMT="(T2,A,T61,A20)") &
2196 0 : "QS| s2 restraint functional form :", "CONSTRAINT"
2197 : CASE DEFAULT
2198 0 : CPABORT("Unknown ddapc restraint")
2199 : END SELECT
2200 : END IF
2201 : END IF
2202 : CALL cp_print_key_finished_output(output_unit, logger, print_section_vals, &
2203 5080 : "DFT_CONTROL_PARAMETERS")
2204 :
2205 5080 : CALL timestop(handle)
2206 :
2207 : END SUBROUTINE write_qs_control
2208 :
2209 : ! **************************************************************************************************
2210 : !> \brief reads the input parameters needed for ddapc.
2211 : !> \param qs_control ...
2212 : !> \param qs_section ...
2213 : !> \param ddapc_restraint_section ...
2214 : !> \author fschiff
2215 : !> \note
2216 : !> either reads DFT%QS%DDAPC_RESTRAINT or PROPERTIES%ET_coupling
2217 : !> if(qs_section is present the DFT part is read, if ddapc_restraint_section
2218 : !> is present ET_COUPLING is read. Avoid having both!!!
2219 : ! **************************************************************************************************
2220 14 : SUBROUTINE read_ddapc_section(qs_control, qs_section, ddapc_restraint_section)
2221 :
2222 : TYPE(qs_control_type), INTENT(INOUT) :: qs_control
2223 : TYPE(section_vals_type), OPTIONAL, POINTER :: qs_section, ddapc_restraint_section
2224 :
2225 : INTEGER :: i, j, jj, k, n_rep
2226 14 : INTEGER, DIMENSION(:), POINTER :: tmplist
2227 14 : REAL(KIND=dp), DIMENSION(:), POINTER :: rtmplist
2228 : TYPE(ddapc_restraint_type), POINTER :: ddapc_restraint_control
2229 : TYPE(section_vals_type), POINTER :: ddapc_section
2230 :
2231 14 : IF (PRESENT(ddapc_restraint_section)) THEN
2232 0 : IF (ASSOCIATED(qs_control%ddapc_restraint_control)) THEN
2233 0 : IF (SIZE(qs_control%ddapc_restraint_control) .GE. 2) &
2234 0 : CPABORT("ET_COUPLING cannot be used in combination with a normal restraint")
2235 : ELSE
2236 0 : ddapc_section => ddapc_restraint_section
2237 0 : ALLOCATE (qs_control%ddapc_restraint_control(1))
2238 : END IF
2239 : END IF
2240 :
2241 14 : IF (PRESENT(qs_section)) THEN
2242 14 : NULLIFY (ddapc_section)
2243 : ddapc_section => section_vals_get_subs_vals(qs_section, &
2244 14 : "DDAPC_RESTRAINT")
2245 : END IF
2246 :
2247 32 : DO i = 1, SIZE(qs_control%ddapc_restraint_control)
2248 :
2249 18 : CALL ddapc_control_create(qs_control%ddapc_restraint_control(i))
2250 18 : ddapc_restraint_control => qs_control%ddapc_restraint_control(i)
2251 :
2252 : CALL section_vals_val_get(ddapc_section, "STRENGTH", i_rep_section=i, &
2253 18 : r_val=ddapc_restraint_control%strength)
2254 : CALL section_vals_val_get(ddapc_section, "TARGET", i_rep_section=i, &
2255 18 : r_val=ddapc_restraint_control%target)
2256 : CALL section_vals_val_get(ddapc_section, "FUNCTIONAL_FORM", i_rep_section=i, &
2257 18 : i_val=ddapc_restraint_control%functional_form)
2258 : CALL section_vals_val_get(ddapc_section, "ATOMS", i_rep_section=i, &
2259 18 : n_rep_val=n_rep)
2260 : CALL section_vals_val_get(ddapc_section, "TYPE_OF_DENSITY", i_rep_section=i, &
2261 18 : i_val=ddapc_restraint_control%density_type)
2262 :
2263 18 : jj = 0
2264 36 : DO k = 1, n_rep
2265 : CALL section_vals_val_get(ddapc_section, "ATOMS", i_rep_section=i, &
2266 18 : i_rep_val=k, i_vals=tmplist)
2267 56 : DO j = 1, SIZE(tmplist)
2268 38 : jj = jj + 1
2269 : END DO
2270 : END DO
2271 18 : IF (jj < 1) CPABORT("Need at least 1 atom to use ddapc constraints")
2272 18 : ddapc_restraint_control%natoms = jj
2273 18 : IF (ASSOCIATED(ddapc_restraint_control%atoms)) &
2274 0 : DEALLOCATE (ddapc_restraint_control%atoms)
2275 54 : ALLOCATE (ddapc_restraint_control%atoms(ddapc_restraint_control%natoms))
2276 18 : jj = 0
2277 36 : DO k = 1, n_rep
2278 : CALL section_vals_val_get(ddapc_section, "ATOMS", i_rep_section=i, &
2279 18 : i_rep_val=k, i_vals=tmplist)
2280 56 : DO j = 1, SIZE(tmplist)
2281 20 : jj = jj + 1
2282 38 : ddapc_restraint_control%atoms(jj) = tmplist(j)
2283 : END DO
2284 : END DO
2285 :
2286 18 : IF (ASSOCIATED(ddapc_restraint_control%coeff)) &
2287 0 : DEALLOCATE (ddapc_restraint_control%coeff)
2288 54 : ALLOCATE (ddapc_restraint_control%coeff(ddapc_restraint_control%natoms))
2289 38 : ddapc_restraint_control%coeff = 1.0_dp
2290 :
2291 : CALL section_vals_val_get(ddapc_section, "COEFF", i_rep_section=i, &
2292 18 : n_rep_val=n_rep)
2293 18 : jj = 0
2294 20 : DO k = 1, n_rep
2295 : CALL section_vals_val_get(ddapc_section, "COEFF", i_rep_section=i, &
2296 2 : i_rep_val=k, r_vals=rtmplist)
2297 22 : DO j = 1, SIZE(rtmplist)
2298 2 : jj = jj + 1
2299 2 : IF (jj > ddapc_restraint_control%natoms) &
2300 0 : CPABORT("Need the same number of coeff as there are atoms ")
2301 4 : ddapc_restraint_control%coeff(jj) = rtmplist(j)
2302 : END DO
2303 : END DO
2304 18 : IF (jj < ddapc_restraint_control%natoms .AND. jj .NE. 0) &
2305 50 : CPABORT("Need no or the same number of coeff as there are atoms.")
2306 : END DO
2307 14 : k = 0
2308 32 : DO i = 1, SIZE(qs_control%ddapc_restraint_control)
2309 18 : IF (qs_control%ddapc_restraint_control(i)%functional_form == &
2310 24 : do_ddapc_constraint) k = k + 1
2311 : END DO
2312 14 : IF (k == 2) CALL cp_abort(__LOCATION__, &
2313 0 : "Only a single constraint possible yet, try to use restraints instead ")
2314 :
2315 14 : END SUBROUTINE read_ddapc_section
2316 :
2317 : ! **************************************************************************************************
2318 : !> \brief ...
2319 : !> \param dft_control ...
2320 : !> \param efield_section ...
2321 : ! **************************************************************************************************
2322 244 : SUBROUTINE read_efield_sections(dft_control, efield_section)
2323 : TYPE(dft_control_type), POINTER :: dft_control
2324 : TYPE(section_vals_type), POINTER :: efield_section
2325 :
2326 : CHARACTER(len=default_path_length) :: file_name
2327 : INTEGER :: i, io, j, n, unit_nr
2328 244 : REAL(KIND=dp), DIMENSION(:), POINTER :: tmp_vals
2329 : TYPE(efield_type), POINTER :: efield
2330 : TYPE(section_vals_type), POINTER :: tmp_section
2331 :
2332 488 : DO i = 1, SIZE(dft_control%efield_fields)
2333 244 : NULLIFY (dft_control%efield_fields(i)%efield)
2334 1220 : ALLOCATE (dft_control%efield_fields(i)%efield)
2335 244 : efield => dft_control%efield_fields(i)%efield
2336 244 : NULLIFY (efield%envelop_i_vars, efield%envelop_r_vars)
2337 : CALL section_vals_val_get(efield_section, "INTENSITY", i_rep_section=i, &
2338 244 : r_val=efield%strength)
2339 :
2340 : CALL section_vals_val_get(efield_section, "POLARISATION", i_rep_section=i, &
2341 244 : r_vals=tmp_vals)
2342 732 : ALLOCATE (efield%polarisation(SIZE(tmp_vals)))
2343 1952 : efield%polarisation = tmp_vals
2344 : CALL section_vals_val_get(efield_section, "PHASE", i_rep_section=i, &
2345 244 : r_val=efield%phase_offset)
2346 : CALL section_vals_val_get(efield_section, "ENVELOP", i_rep_section=i, &
2347 244 : i_val=efield%envelop_id)
2348 : CALL section_vals_val_get(efield_section, "WAVELENGTH", i_rep_section=i, &
2349 244 : r_val=efield%wavelength)
2350 : CALL section_vals_val_get(efield_section, "VEC_POT_INITIAL", i_rep_section=i, &
2351 244 : r_vals=tmp_vals)
2352 1952 : efield%vec_pot_initial = tmp_vals
2353 :
2354 488 : IF (efield%envelop_id == constant_env) THEN
2355 238 : ALLOCATE (efield%envelop_i_vars(2))
2356 238 : tmp_section => section_vals_get_subs_vals(efield_section, "CONSTANT_ENV", i_rep_section=i)
2357 : CALL section_vals_val_get(tmp_section, "START_STEP", &
2358 238 : i_val=efield%envelop_i_vars(1))
2359 : CALL section_vals_val_get(tmp_section, "END_STEP", &
2360 238 : i_val=efield%envelop_i_vars(2))
2361 6 : ELSE IF (efield%envelop_id == gaussian_env) THEN
2362 2 : ALLOCATE (efield%envelop_r_vars(2))
2363 2 : tmp_section => section_vals_get_subs_vals(efield_section, "GAUSSIAN_ENV", i_rep_section=i)
2364 : CALL section_vals_val_get(tmp_section, "T0", &
2365 2 : r_val=efield%envelop_r_vars(1))
2366 : CALL section_vals_val_get(tmp_section, "SIGMA", &
2367 2 : r_val=efield%envelop_r_vars(2))
2368 4 : ELSE IF (efield%envelop_id == ramp_env) THEN
2369 2 : ALLOCATE (efield%envelop_i_vars(4))
2370 2 : tmp_section => section_vals_get_subs_vals(efield_section, "RAMP_ENV", i_rep_section=i)
2371 : CALL section_vals_val_get(tmp_section, "START_STEP_IN", &
2372 2 : i_val=efield%envelop_i_vars(1))
2373 : CALL section_vals_val_get(tmp_section, "END_STEP_IN", &
2374 2 : i_val=efield%envelop_i_vars(2))
2375 : CALL section_vals_val_get(tmp_section, "START_STEP_OUT", &
2376 2 : i_val=efield%envelop_i_vars(3))
2377 : CALL section_vals_val_get(tmp_section, "END_STEP_OUT", &
2378 2 : i_val=efield%envelop_i_vars(4))
2379 2 : ELSE IF (efield%envelop_id == custom_env) THEN
2380 2 : tmp_section => section_vals_get_subs_vals(efield_section, "CUSTOM_ENV", i_rep_section=i)
2381 2 : CALL section_vals_val_get(tmp_section, "EFIELD_FILE_NAME", c_val=file_name)
2382 2 : CALL open_file(file_name=TRIM(file_name), file_action="READ", file_status="OLD", unit_number=unit_nr)
2383 : !Determine the number of lines in file
2384 2 : n = 0
2385 10 : DO WHILE (.TRUE.)
2386 12 : READ (unit_nr, *, iostat=io)
2387 12 : IF (io /= 0) EXIT
2388 10 : n = n + 1
2389 : END DO
2390 2 : REWIND (unit_nr)
2391 6 : ALLOCATE (efield%envelop_r_vars(n + 1))
2392 : !Store the timestep of the list in the first entry of the r_vars
2393 2 : CALL section_vals_val_get(tmp_section, "TIMESTEP", r_val=efield%envelop_r_vars(1))
2394 : !Read the file
2395 12 : DO j = 2, n + 1
2396 10 : READ (unit_nr, *) efield%envelop_r_vars(j)
2397 12 : efield%envelop_r_vars(j) = cp_unit_to_cp2k(efield%envelop_r_vars(j), "volt/m")
2398 : END DO
2399 2 : CALL close_file(unit_nr)
2400 : END IF
2401 : END DO
2402 244 : END SUBROUTINE read_efield_sections
2403 :
2404 : ! **************************************************************************************************
2405 : !> \brief reads the input parameters needed real time propagation
2406 : !> \param dft_control ...
2407 : !> \param rtp_section ...
2408 : !> \author fschiff
2409 : ! **************************************************************************************************
2410 472 : SUBROUTINE read_rtp_section(dft_control, rtp_section)
2411 :
2412 : TYPE(dft_control_type), INTENT(INOUT) :: dft_control
2413 : TYPE(section_vals_type), POINTER :: rtp_section
2414 :
2415 236 : INTEGER, DIMENSION(:), POINTER :: tmp
2416 : LOGICAL :: is_present
2417 : TYPE(section_vals_type), POINTER :: proj_mo_section
2418 :
2419 2832 : ALLOCATE (dft_control%rtp_control)
2420 : CALL section_vals_val_get(rtp_section, "MAX_ITER", &
2421 236 : i_val=dft_control%rtp_control%max_iter)
2422 : CALL section_vals_val_get(rtp_section, "MAT_EXP", &
2423 236 : i_val=dft_control%rtp_control%mat_exp)
2424 : CALL section_vals_val_get(rtp_section, "ASPC_ORDER", &
2425 236 : i_val=dft_control%rtp_control%aspc_order)
2426 : CALL section_vals_val_get(rtp_section, "EXP_ACCURACY", &
2427 236 : r_val=dft_control%rtp_control%eps_exp)
2428 : CALL section_vals_val_get(rtp_section, "PROPAGATOR", &
2429 236 : i_val=dft_control%rtp_control%propagator)
2430 : CALL section_vals_val_get(rtp_section, "EPS_ITER", &
2431 236 : r_val=dft_control%rtp_control%eps_ener)
2432 : CALL section_vals_val_get(rtp_section, "INITIAL_WFN", &
2433 236 : i_val=dft_control%rtp_control%initial_wfn)
2434 : CALL section_vals_val_get(rtp_section, "HFX_BALANCE_IN_CORE", &
2435 236 : l_val=dft_control%rtp_control%hfx_redistribute)
2436 : CALL section_vals_val_get(rtp_section, "APPLY_WFN_MIX_INIT_RESTART", &
2437 236 : l_val=dft_control%rtp_control%apply_wfn_mix_init_restart)
2438 : CALL section_vals_val_get(rtp_section, "APPLY_DELTA_PULSE", &
2439 236 : l_val=dft_control%rtp_control%apply_delta_pulse)
2440 : CALL section_vals_val_get(rtp_section, "APPLY_DELTA_PULSE_MAG", &
2441 236 : l_val=dft_control%rtp_control%apply_delta_pulse_mag)
2442 : CALL section_vals_val_get(rtp_section, "VELOCITY_GAUGE", &
2443 236 : l_val=dft_control%rtp_control%velocity_gauge)
2444 : CALL section_vals_val_get(rtp_section, "VG_COM_NL", &
2445 236 : l_val=dft_control%rtp_control%nl_gauge_transform)
2446 : CALL section_vals_val_get(rtp_section, "PERIODIC", &
2447 236 : l_val=dft_control%rtp_control%periodic)
2448 : CALL section_vals_val_get(rtp_section, "DENSITY_PROPAGATION", &
2449 236 : l_val=dft_control%rtp_control%linear_scaling)
2450 : CALL section_vals_val_get(rtp_section, "MCWEENY_MAX_ITER", &
2451 236 : i_val=dft_control%rtp_control%mcweeny_max_iter)
2452 : CALL section_vals_val_get(rtp_section, "ACCURACY_REFINEMENT", &
2453 236 : i_val=dft_control%rtp_control%acc_ref)
2454 : CALL section_vals_val_get(rtp_section, "MCWEENY_EPS", &
2455 236 : r_val=dft_control%rtp_control%mcweeny_eps)
2456 : CALL section_vals_val_get(rtp_section, "DELTA_PULSE_SCALE", &
2457 236 : r_val=dft_control%rtp_control%delta_pulse_scale)
2458 : CALL section_vals_val_get(rtp_section, "DELTA_PULSE_DIRECTION", &
2459 236 : i_vals=tmp)
2460 944 : dft_control%rtp_control%delta_pulse_direction = tmp
2461 : CALL section_vals_val_get(rtp_section, "SC_CHECK_START", &
2462 236 : i_val=dft_control%rtp_control%sc_check_start)
2463 236 : proj_mo_section => section_vals_get_subs_vals(rtp_section, "PRINT%PROJECTION_MO")
2464 236 : CALL section_vals_get(proj_mo_section, explicit=is_present)
2465 236 : IF (is_present) THEN
2466 4 : IF (dft_control%rtp_control%linear_scaling) &
2467 : CALL cp_abort(__LOCATION__, &
2468 : "You have defined a time dependent projection of mos, but "// &
2469 : "only the density matrix is propagated (DENSITY_PROPAGATION "// &
2470 : ".TRUE.). Please either use MO-based real time DFT or do not "// &
2471 0 : "define any PRINT%PROJECTION_MO section")
2472 4 : dft_control%rtp_control%is_proj_mo = .TRUE.
2473 : ELSE
2474 232 : dft_control%rtp_control%is_proj_mo = .FALSE.
2475 : END IF
2476 :
2477 236 : END SUBROUTINE read_rtp_section
2478 :
2479 : ! **************************************************************************************************
2480 : !> \brief Parses the BLOCK_LIST keywords from the ADMM section
2481 : !> \param admm_control ...
2482 : !> \param dft_section ...
2483 : ! **************************************************************************************************
2484 442 : SUBROUTINE read_admm_block_list(admm_control, dft_section)
2485 : TYPE(admm_control_type), POINTER :: admm_control
2486 : TYPE(section_vals_type), POINTER :: dft_section
2487 :
2488 : INTEGER :: irep, list_size, n_rep
2489 442 : INTEGER, DIMENSION(:), POINTER :: tmplist
2490 :
2491 442 : NULLIFY (tmplist)
2492 :
2493 : CALL section_vals_val_get(dft_section, "AUXILIARY_DENSITY_MATRIX_METHOD%BLOCK_LIST", &
2494 442 : n_rep_val=n_rep)
2495 :
2496 938 : ALLOCATE (admm_control%blocks(n_rep))
2497 :
2498 478 : DO irep = 1, n_rep
2499 : CALL section_vals_val_get(dft_section, "AUXILIARY_DENSITY_MATRIX_METHOD%BLOCK_LIST", &
2500 36 : i_rep_val=irep, i_vals=tmplist)
2501 36 : list_size = SIZE(tmplist)
2502 108 : ALLOCATE (admm_control%blocks(irep)%list(list_size))
2503 650 : admm_control%blocks(irep)%list(:) = tmplist(:)
2504 : END DO
2505 :
2506 442 : END SUBROUTINE read_admm_block_list
2507 :
2508 : END MODULE cp_control_utils
|