35#include "./base/base_uses.f90"
40 CHARACTER(len=*),
PARAMETER,
PRIVATE :: moduleN =
'input_cp2k_se'
57 cpassert(.NOT.
ASSOCIATED(section))
59 description=
"Parameters needed to set up the Semi-empirical methods", &
60 n_keywords=8, n_subsections=0, repeats=.false.)
62 CALL keyword_create(keyword, __location__, name=
"ORTHOGONAL_BASIS", &
63 description=
"Assume orthogonal basis set. This flag is overwritten by "// &
64 "methods with fixed orthogonal/non-orthogonal basis set.", &
65 usage=
"ORTHOGONAL_BASIS", default_l_val=.false., lone_keyword_l_val=.true.)
70 description=
"Provides the order of the Slater orbital expansion of Gaussian-Type Orbitals.", &
71 usage=
"STO_NG", default_i_val=6)
75 CALL keyword_create(keyword, __location__, name=
"ANALYTICAL_GRADIENTS", &
76 description=
"Nuclear Gradients are computed analytically or numerically", &
77 usage=
"ANALYTICAL_GRADIENTS", default_l_val=.true., lone_keyword_l_val=.true.)
82 description=
"Step size in finite difference force calculation", &
83 usage=
"DELTA {real} ", default_r_val=1.e-6_dp)
87 CALL keyword_create(keyword, __location__, name=
"INTEGRAL_SCREENING", &
88 description=
"Specifies the functional form for the ", &
89 usage=
"INTEGRAL_SCREENING (KDSO|KDSO-D|SLATER)", &
90 enum_c_vals=
s2a(
"KDSO",
"KDSO-D",
"SLATER"), &
92 enum_desc=
s2a(
"Uses the standard NDDO Klopman-Dewar-Sabelli-Ohno equation "// &
93 "for the screening of the Coulomb interactions.", &
94 "Uses a modified Klopman-Dewar-Sabelli-Ohno equation, dumping the screening "// &
95 "parameter for the Coulomb interactions.", &
96 "Uses an exponential Slater-type function for modelling the Coulomb interactions."), &
102 description=
"Specifies the type of treatment for the electrostatic long-range part "// &
103 "in semi-empirical calculations.", &
104 usage=
"PERIODIC (NONE|EWALD|EWALD_R3|EWALD_GKS)", &
105 enum_c_vals=
s2a(
"NONE",
"EWALD",
"EWALD_R3",
"EWALD_GKS"), &
107 enum_desc=
s2a(
"The long-range part is not explicitly treaten. The interaction "// &
108 "depends uniquely on the Cutoffs used for the calculation.", &
109 "Enables the usage of Multipoles Ewald summation schemes. The short-range part "// &
110 "is tapered according the value of RC_COULOMB.", &
111 "Enables the usage of Multipoles Ewald summation schemes together with a long-range "// &
112 "treatment for the 1/R^3 term, deriving from the short-range component. This option "// &
113 "is active only for K-DSO integral screening type.", &
114 "Use Ewald directly in Coulomb integral evaluation, works only with Slater screening"), &
119 CALL keyword_create(keyword, __location__, name=
"FORCE_KDSO-D_EXCHANGE", &
120 description=
"This keywords forces the usage of the KDSO-D integral screening "// &
121 "for the Exchange integrals (default is to apply the screening only to the "// &
122 "Coulomb integrals.", default_l_val=.false., lone_keyword_l_val=.true.)
127 description=
"Use dispersion correction", &
128 lone_keyword_l_val=.true., &
129 usage=
"DISPERSION", default_l_val=.false.)
133 CALL keyword_create(keyword, __location__, name=
"DISPERSION_PARAMETER_FILE", &
134 description=
"Specify file that contains the atomic dispersion parameters", &
135 usage=
"DISPERSION_PARAMETER_FILE filename", &
136 n_var=1, type_of_var=
char_t, default_c_val=
"")
140 CALL keyword_create(keyword, __location__, name=
"DISPERSION_RADIUS", &
141 description=
"Define radius of dispersion interaction", &
142 usage=
"DISPERSION_RADIUS", default_r_val=15._dp)
146 CALL keyword_create(keyword, __location__, name=
"COORDINATION_CUTOFF", &
147 description=
"Define cutoff for coordination number calculation", &
148 usage=
"COORDINATION_CUTOFF", default_r_val=1.e-6_dp)
153 description=
"Scaling parameters (s6,sr6,s8) for the D3 dispersion method,", &
154 usage=
"D3_SCALING 1.0 1.0 1.0", n_var=3, default_r_vals=(/0.0_dp, 0.0_dp, 0.0_dp/))
159 CALL create_coulomb_section(subsection)
163 CALL create_exchange_section(subsection)
167 CALL create_screening_section(subsection)
171 CALL create_lr_corr_section(subsection)
179 CALL create_se_memory_section(subsection)
183 CALL create_se_print_section(subsection)
187 CALL create_se_ga_section(subsection)
199 SUBROUTINE create_lr_corr_section(section)
204 cpassert(.NOT.
ASSOCIATED(section))
205 CALL section_create(section, __location__, name=
"LR_CORRECTION", &
206 description=
"Setup parameters for the evaluation of the long-range correction term in SE "// &
207 "calculations.", n_keywords=0, n_subsections=1, repeats=.false.)
211 description=
"Atomic Cutoff Radius Cutoff for the evaluation of the long-ranbe correction integrals. ", &
212 usage=
"CUTOFF {real} ", unit_str=
"angstrom", &
218 description=
"Atomic Cutoff Radius Cutoff for Tapering the long-range correction integrals. "// &
219 "If not specified it assumes the same value specified for the CUTOFF.", &
220 usage=
"RC_TAPER {real} ", unit_str=
"angstrom", type_of_var=
real_t)
225 description=
"Range of cutoff switch function (tapering): 0.5*(1-TANH((r-r0)/RC_RANGE)), "// &
226 "where r0=2.0*RC_TAPER-20.0*RC_RANGE.", &
227 usage=
"RC_RANGE {real} ", unit_str=
"angstrom", default_r_val=0.0_dp)
231 END SUBROUTINE create_lr_corr_section
239 SUBROUTINE create_coulomb_section(section)
244 cpassert(.NOT.
ASSOCIATED(section))
246 description=
"Setup parameters for the evaluation of the COULOMB term in SE "// &
247 "calculations.", n_keywords=0, n_subsections=1, repeats=.false.)
251 keyword, __location__, name=
"CUTOFF", &
252 description=
"Atomic Cutoff Radius Cutoff for the evaluation of the Coulomb integrals. "// &
253 "For non-periodic calculation the default value is exactly the full cell dimension, in order "// &
254 "to evaluate all pair interactions. Instead, for periodic calculations the default numerical value is used.", &
255 usage=
"CUTOFF {real} ", unit_str=
"angstrom", &
261 description=
"Atomic Cutoff Radius Cutoff for Tapering Coulomb integrals. "// &
262 "If not specified it assumes the same value specified for the CUTOFF.", &
263 usage=
"RC_TAPER {real} ", unit_str=
"angstrom", type_of_var=
real_t)
268 description=
"Range of cutoff switch function (tapering): 0.5*(1-TANH((r-r0)/RC_RANGE)), "// &
269 "where r0=2.0*RC_TAPER-20.0*RC_RANGE.", &
270 usage=
"RC_RANGE {real} ", unit_str=
"angstrom", default_r_val=0.0_dp)
274 END SUBROUTINE create_coulomb_section
282 SUBROUTINE create_exchange_section(section)
287 cpassert(.NOT.
ASSOCIATED(section))
289 description=
"Setup parameters for the evaluation of the EXCHANGE and "// &
290 "core Hamiltonian terms in SE calculations.", n_keywords=0, n_subsections=1, &
295 description=
"Atomic Cutoff Radius Cutoff for the evaluation of the Exchange integrals. "// &
296 "For non-periodic calculation the default value is exactly the full cell dimension, in order "// &
297 "to evaluate all pair interactions. Instead, for periodic calculations the default is the "// &
298 "minimum value between 1/4 of the cell dimension and the value specified in input (either"// &
299 " explicitly defined or the default numerical value).", &
300 usage=
"CUTOFF {real} ", unit_str=
"angstrom", &
306 description=
"Atomic Cutoff Radius Cutoff for Tapering Exchange integrals. "// &
307 "If not specified it assumes the same value specified for the CUTOFF.", &
308 usage=
"RC_TAPER {real} ", unit_str=
"angstrom", type_of_var=
real_t)
313 description=
"Range of cutoff switch function (tapering): 0.5*(1-TANH((r-r0)/RC_RANGE)), "// &
314 "where r0=2.0*RC_TAPER-20.0*RC_RANGE.", &
315 usage=
"RC_RANGE {real} ", unit_str=
"angstrom", default_r_val=0.0_dp)
319 END SUBROUTINE create_exchange_section
327 SUBROUTINE create_screening_section(section)
332 cpassert(.NOT.
ASSOCIATED(section))
334 description=
"Setup parameters for the tapering of the Coulomb/Exchange Screening in "// &
335 "KDSO-D integral scheme,", n_keywords=0, n_subsections=1, repeats=.false.)
339 description=
"Atomic Cutoff Radius Cutoff for Tapering the screening term. ", &
340 usage=
"RC_TAPER {real} ", unit_str=
"angstrom", &
346 description=
"Range of cutoff switch function (tapering): 0.5*(1-TANH((r-r0)/RC_RANGE)), "// &
347 "where r0=2*RC_TAPER-20*RC_RANGE.", &
348 usage=
"RC_RANGE {real} ", unit_str=
"angstrom", default_r_val=0.0_dp)
352 END SUBROUTINE create_screening_section
359 SUBROUTINE create_se_print_section(section)
364 cpassert(.NOT.
ASSOCIATED(section))
366 description=
"Section of possible print options in SE code.", &
367 n_keywords=0, n_subsections=1, repeats=.false.)
371 description=
"Activates the printing of the neighbor lists used"// &
372 " for the periodic SE calculations.", &
378 description=
"Activates the printing of the subcells used for the "// &
379 "generation of neighbor lists for periodic SE.", &
385 description=
"Activates the printing of the information for "// &
386 "Ewald multipole summation in periodic SE.", &
391 END SUBROUTINE create_se_print_section
398 SUBROUTINE create_se_ga_section(section)
403 cpassert(.NOT.
ASSOCIATED(section))
405 description=
"Sets up memory parameters for the storage of the integrals", &
406 n_keywords=1, n_subsections=0, repeats=.false.)
410 keyword, __location__, &
412 description=
"Defines the number of linked cells for the neighbor list. "// &
413 "Default value is number of processors", &
418 END SUBROUTINE create_se_ga_section
425 SUBROUTINE create_se_memory_section(section)
430 cpassert(.NOT.
ASSOCIATED(section))
432 description=
"Sets up memory parameters for the storage of the integrals", &
433 n_keywords=1, n_subsections=0, repeats=.false.)
436 keyword, __location__, &
437 name=
"EPS_STORAGE", &
438 description=
"Storage threshold for compression is EPS_STORAGE", &
439 usage=
"EPS_STORAGE 1.0E-10", &
440 default_r_val=1.0e-10_dp)
445 keyword, __location__, &
447 description=
"Defines the maximum amount of memory [MB] used to store precomputed "// &
448 "(possibly compressed) two-electron two-center integrals", &
449 usage=
"MAX_MEMORY 256", &
455 description=
"Enables the compression of the integrals in memory.", &
456 usage=
"COMPRESS <LOGICAL>", &
457 default_l_val=.false., lone_keyword_l_val=.true.)
461 END SUBROUTINE create_se_memory_section
routines to handle the output, The idea is to remove the decision of wheter to output and what to out...
integer, parameter, public high_print_level
subroutine, public cp_print_key_section_create(print_key_section, location, name, description, print_level, each_iter_names, each_iter_values, add_last, filename, common_iter_levels, citations, unit_str)
creates a print_key section
real(kind=dp) function, public cp_unit_to_cp2k(value, unit_str, defaults, power)
converts to the internal cp2k units to the given unit
Defines the basic variable types.
integer, parameter, public dp
Utilities for string manipulations.