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