(git:b279b6b)
input_cp2k_global.F
Go to the documentation of this file.
1 !--------------------------------------------------------------------------------------------------!
2 ! CP2K: A general program to perform molecular dynamics simulations !
3 ! Copyright 2000-2024 CP2K developers group <https://cp2k.org> !
4 ! !
5 ! SPDX-License-Identifier: GPL-2.0-or-later !
6 !--------------------------------------------------------------------------------------------------!
7 
8 ! **************************************************************************************************
9 !> \brief builds the global input section for cp2k
10 !> \par History
11 !> 06.2004 created [fawzi]
12 !> 03.2014 moved to separate file [Ole Schuett]
13 !> 10.2016 update seed input [Matthias Krack]
14 !> \author fawzi
15 ! **************************************************************************************************
17  USE bibliography, ONLY: ceriotti2014,&
18  frigo2005,&
20  USE cp_blacs_env, ONLY: blacs_grid_col,&
39  USE grid_api, ONLY: grid_backend_auto,&
45  USE input_constants, ONLY: &
55  keyword_type
60  section_type
61  USE input_val_types, ONLY: char_t,&
62  integer_t,&
63  logical_t
64  USE kinds, ONLY: dp
65  USE string_utilities, ONLY: s2a
67 #include "./base/base_uses.f90"
68 
69  IMPLICIT NONE
70  PRIVATE
71 
72  LOGICAL, PRIVATE, PARAMETER :: debug_this_module = .true.
73  CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'input_cp2k_global'
74 
75  PUBLIC :: create_global_section
76 
77 CONTAINS
78 
79 ! **************************************************************************************************
80 !> \brief section to hold global settings for the whole program
81 !> \param section the section to be created
82 !> \author fawzi
83 ! **************************************************************************************************
84  SUBROUTINE create_global_section(section)
85  TYPE(section_type), POINTER :: section
86 
87  INTEGER :: default_dgemm
88  TYPE(keyword_type), POINTER :: keyword
89  TYPE(section_type), POINTER :: print_key, sub_section
90 
91  NULLIFY (print_key)
92  cpassert(.NOT. ASSOCIATED(section))
93  CALL section_create(section, __location__, name="GLOBAL", &
94  description="Section with general information regarding which kind "// &
95  "of simulation to perform an parameters for the whole PROGRAM", &
96  n_keywords=7, n_subsections=0, repeats=.false.)
97 
98  NULLIFY (keyword)
99  CALL keyword_create(keyword, __location__, name="BLACS_GRID", &
100  description="how to distribute the processors on the 2d grid needed "// &
101  "by BLACS (and thus SCALAPACK)", usage="BLACS_GRID SQUARE", &
102  default_i_val=blacs_grid_square, enum_c_vals=s2a("SQUARE", "ROW", "COLUMN"), &
103  enum_desc=s2a("Distribution by matrix blocks", "Distribution by matrix rows", &
104  "Distribution by matrix columns"), &
106  CALL section_add_keyword(section, keyword)
107  CALL keyword_release(keyword)
108 
109  CALL keyword_create(keyword, __location__, name="BLACS_REPEATABLE", &
110  description="Use a topology for BLACS collectives that is guaranteed to be repeatable "// &
111  "on homogeneous architectures", &
112  usage="BLACS_REPEATABLE", &
113  default_l_val=.false., lone_keyword_l_val=.true.)
114  CALL section_add_keyword(section, keyword)
115  CALL keyword_release(keyword)
116 
117  CALL keyword_create(keyword, __location__, name="PREFERRED_DIAG_LIBRARY", &
118  description="Specifies the diagonalization library to be used. If not available, "// &
119  "the ScaLAPACK library is used", &
120  usage="PREFERRED_DIAG_LIBRARY ELPA", &
121  enum_i_vals=(/fm_diag_type_elpa, &
125  fm_diag_type_dlaf/), &
126  enum_c_vals=s2a("ELPA", "ScaLAPACK", "SL", "CUSOLVER", "DLAF"), &
127  enum_desc=s2a("ELPA library", &
128  "ScaLAPACK library", &
129  "ScaLAPACK library (shorthand)", &
130  "cuSOLVER (CUDA GPU library)", &
131  "DLA-Future (CUDA/HIP GPU library)"), &
132  default_i_val=fm_diag_type_default)
133  CALL section_add_keyword(section, keyword)
134  CALL keyword_release(keyword)
135 
136 #if defined(__SPLA) && defined(__OFFLOAD_GEMM)
137  default_dgemm = do_dgemm_spla
138 #else
139  default_dgemm = do_dgemm_blas
140 #endif
141  CALL keyword_create(keyword, __location__, name="PREFERRED_DGEMM_LIBRARY", &
142  description="Specifies the DGEMM library to be used. If not available, "// &
143  "the BLAS routine is used. This keyword affects some DGEMM calls in the WFC code and turns on their "// &
144  "acceleration with SpLA. This keyword affects only local DGEMM calls, not the calls to PDGEMM "// &
145  "(see keyword FM%TYPE_OF_MATRIX_MULTIPLICATION).", &
146  usage="PREFERRED_DGEMM_LIBRARY SPLA", &
147  default_i_val=default_dgemm, &
148  enum_i_vals=(/do_dgemm_spla, do_dgemm_blas/), &
149  enum_c_vals=s2a("SPLA", "BLAS"), &
150  enum_desc=s2a("SPLA library", "BLAS library"))
151  CALL section_add_keyword(section, keyword)
152  CALL keyword_release(keyword)
153 
154  CALL keyword_create(keyword, __location__, name="EPS_CHECK_DIAG", &
155  description="Check that the orthonormality of the eigenvectors after a diagonalization "// &
156  "fulfills the specified numerical accuracy. A negative threshold value disables the check.", &
157  usage="EPS_CHECK_DIAG 1.0E-14", &
158  default_r_val=-1.0_dp)
159  CALL section_add_keyword(section, keyword)
160  CALL keyword_release(keyword)
161 
162  CALL keyword_create(keyword, __location__, name="ELPA_KERNEL", &
163  description="Specifies the kernel to be used when ELPA is in use", &
164  default_i_val=elpa_kernel_ids(1), &
165  enum_i_vals=elpa_kernel_ids, &
166  enum_c_vals=elpa_kernel_names, &
167  enum_desc=elpa_kernel_descriptions)
168  CALL section_add_keyword(section, keyword)
169  CALL keyword_release(keyword)
170 
171  CALL keyword_create(keyword, __location__, name="ELPA_NEIGVEC_MIN", &
172  description="Minimum number of eigenvectors for the use of the eigensolver from "// &
173  "the ELPA library. The eigensolver from the ScaLAPACK library is used as fallback "// &
174  "for all smaller cases", &
175  usage="ELPA_NEIGVEC_MIN 32", &
176  default_i_val=16)
177  CALL section_add_keyword(section, keyword)
178  CALL keyword_release(keyword)
179 
180  CALL keyword_create(keyword, __location__, name="ELPA_QR", &
181  description="For ELPA, enable a blocked QR step when reducing the input matrix "// &
182  "to banded form in preparation for the actual diagonalization step. "// &
183  "See implementation paper for more details. Requires ELPA version 201505 or newer, "// &
184  "automatically deactivated otherwise. If true, QR is activated only when the "// &
185  "the size of the diagonalized matrix is suitable. Print key PRINT_ELPA is "// &
186  "useful in determining which matrices are suitable for QR. Might accelerate the "// &
187  "diagonalization of suitable matrices.", &
188  usage="ELPA_QR", &
189  default_l_val=.false., lone_keyword_l_val=.true.)
190  CALL section_add_keyword(section, keyword)
191  CALL keyword_release(keyword)
192 
193  CALL keyword_create(keyword, __location__, name="ELPA_QR_UNSAFE", &
194  description="For ELPA, disable block size limitations when used together with ELPA_QR. "// &
195  "Keyword relevant only with ELPA versions 201605 or newer. Use keyword with caution, "// &
196  "as it might result in wrong eigenvalues with some matrix orders/block sizes "// &
197  "when the number of MPI processes is varied. If the print key PRINT_ELPA is "// &
198  "active the validity of the eigenvalues is checked against values calculated without "// &
199  "ELPA QR.", &
200  usage="ELPA_QR", &
201  default_l_val=.false., lone_keyword_l_val=.true.)
202  CALL section_add_keyword(section, keyword)
203  CALL keyword_release(keyword)
204 
205  CALL cp_print_key_section_create(print_key, __location__, "PRINT_ELPA", &
206  description="Controls the printing of ELPA diagonalization information. "// &
207  "Useful for testing purposes, especially together with keyword ELPA_QR.", &
208  filename="__STD_OUT__")
209  CALL section_add_subsection(section, print_key)
210  CALL section_release(print_key)
211 
212  CALL keyword_create(keyword, __location__, name="DLAF_NEIGVEC_MIN", &
213  description="Minimum number of eigenvectors for the use of the eigensolver from "// &
214  "the DLA-Future library. The eigensolver from the ScaLAPACK library is used as fallback "// &
215  "for all smaller cases", &
216  usage="DLAF_NEIGVEC_MIN 512", &
217  default_i_val=1024)
218  CALL section_add_keyword(section, keyword)
219  CALL keyword_release(keyword)
220 
221  CALL keyword_create( &
222  keyword, __location__, name="PREFERRED_FFT_LIBRARY", &
223  description="Specifies the FFT library which should be preferred. "// &
224  "If it is not available, use FFTW3 if this is linked in, if FFTW3 is not available use FFTSG. "// &
225  "Improved performance with FFTW3 can be obtained specifying a proper value for FFTW_PLAN_TYPE. "// &
226  "Contrary to earlier CP2K versions, all libraries will result in the same grids, "// &
227  "i.e. the subset of grids which all FFT libraries can transform. "// &
228  "See EXTENDED_FFT_LENGTHS if larger FFTs or grids that more precisely match a given cutoff are needed, "// &
229  "or older results need to be reproduced. "// &
230  "FFTW3 is often (close to) optimal, and well tested with CP2K.", &
231  usage="PREFERRED_FFT_LIBRARY FFTW3", &
232  citations=(/frigo2005/), &
233  default_i_val=do_fft_fftw3, &
234  enum_i_vals=(/do_fft_sg, do_fft_fftw3, do_fft_fftw3/), &
235  enum_c_vals=s2a("FFTSG", "FFTW3", "FFTW"), &
236  enum_desc=s2a("Stefan Goedecker's FFT (FFTSG), always available, "// &
237  "will be used in case a FFT library is specified and not available.", &
238  "a fast portable FFT library. Recommended. "// &
239  "See also the FFTW_PLAN_TYPE, and FFTW_WISDOM_FILE_NAME keywords.", &
240  "Same as FFTW3 (for compatibility with CP2K 2.3)"))
241  CALL section_add_keyword(section, keyword)
242  CALL keyword_release(keyword)
243 
244  CALL keyword_create(keyword, __location__, name="FFTW_WISDOM_FILE_NAME", &
245  description="The name of the file that contains wisdom (pre-planned FFTs) for use with FFTW3. "// &
246  "Using wisdom can significantly speed up the FFTs (see the FFTW homepage for details). "// &
247  "Note that wisdom is not transferable between different computer (architectures). "// &
248  "Wisdom can be generated using the fftw-wisdom tool that is part of the fftw installation. "// &
249  "cp2k/tools/cp2k-wisdom is a script that contains some additional info, and can help "// &
250  "to generate a useful default for /etc/fftw/wisdom or particular values for a given simulation.", &
251  usage="FFTW_WISDOM_FILE_NAME wisdom.dat", default_lc_val="/etc/fftw/wisdom")
252  CALL section_add_keyword(section, keyword)
253  CALL keyword_release(keyword)
254 
255  CALL keyword_create(keyword, __location__, name="FFTW_PLAN_TYPE", &
256  description="FFTW can have improved performance if it is allowed to plan with "// &
257  "explicit measurements which strategy is best for a given FFT. "// &
258  "While a plan based on measurements is generally faster, "// &
259  "differences in machine load will lead to different plans for the same input file, "// &
260  "and thus numerics for the FFTs will be slightly different from run to run. "// &
261  "PATIENT planning is recommended for long ab initio MD runs.", &
262  usage="FFTW_PLAN_TYPE PATIENT", &
263  citations=(/frigo2005/), &
264  default_i_val=fftw_plan_estimate, &
266  enum_c_vals=s2a("ESTIMATE", &
267  "MEASURE", &
268  "PATIENT", &
269  "EXHAUSTIVE"), &
270  enum_desc=s2a("Quick estimate, no runtime measurements.", &
271  "Quick measurement, somewhat faster FFTs.", &
272  "Measurements trying a wider range of possibilities.", &
273  "Measurements trying all possibilities - use with caution."))
274  CALL section_add_keyword(section, keyword)
275  CALL keyword_release(keyword)
276 
277  CALL keyword_create(keyword, __location__, name="EXTENDED_FFT_LENGTHS", &
278  description="Use fft library specific values for the allows number of points in FFTs. "// &
279  "The default is to use the internal FFT lengths. For external fft libraries this may "// &
280  "create an error at the external library level, because the length provided by cp2k is "// &
281  "not supported by the external library. In this case switch on this keyword "// &
282  "to obtain, with certain fft libraries, lengths matching the external fft library lengths, or "// &
283  "larger allowed grids, or grids that more precisely match a given cutoff. "// &
284  "IMPORTANT NOTE: in this case, the actual grids used in CP2K depends on the FFT library. "// &
285  "A change of FFT library must therefore be considered equivalent to a change of basis, "// &
286  "which implies a change of total energy. ", &
287  usage="EXTENDED_FFT_LENGTHS", &
288  default_l_val=.false., lone_keyword_l_val=.true.)
289  CALL section_add_keyword(section, keyword)
290  CALL keyword_release(keyword)
291 
292  CALL keyword_create(keyword, __location__, name="FFT_POOL_SCRATCH_LIMIT", &
293  description="Limits the memory usage of the FFT scratch pool, potentially reducing efficiency a bit", &
294  usage="FFT_POOL_SCRATCH_LIMIT {INTEGER}", default_i_val=15)
295  CALL section_add_keyword(section, keyword)
296  CALL keyword_release(keyword)
297 
298  CALL keyword_create(keyword, __location__, name="ALLTOALL_SGL", &
299  description="All-to-all communication (FFT) should use single precision", &
300  usage="ALLTOALL_SGL YES", &
301  default_l_val=.false., lone_keyword_l_val=.true.)
302  CALL section_add_keyword(section, keyword)
303  CALL keyword_release(keyword)
304 
305  CALL keyword_create(keyword, __location__, name="PRINT_LEVEL", &
306  variants=(/"IOLEVEL"/), &
307  description="How much output is written out.", &
308  usage="PRINT_LEVEL HIGH", &
309  default_i_val=medium_print_level, enum_c_vals= &
310  s2a("SILENT", "LOW", "MEDIUM", "HIGH", "DEBUG"), &
311  enum_desc=s2a("Almost no output", &
312  "Little output", "Quite some output", "Lots of output", &
313  "Everything is written out, useful for debugging purposes only"), &
316  CALL section_add_keyword(section, keyword)
317  CALL keyword_release(keyword)
318 
319  CALL keyword_create( &
320  keyword, __location__, name="PROGRAM_NAME", &
321  variants=(/"PROGRAM"/), &
322  description="Which program should be run", &
323  usage="PROGRAM_NAME {STRING}", &
324  enum_c_vals=s2a("ATOM", "FARMING", "TEST", "CP2K", "OPTIMIZE_INPUT", "OPTIMIZE_BASIS", "TMC", "MC_ANALYSIS", "SWARM"), &
325  enum_desc=s2a("Runs single atom calculations", &
326  "Runs N independent jobs in a single run", &
327  "Do some benchmarking and testing", &
328  "Runs one of the CP2K package", &
329  "A tool to optimize parameters in a CP2K input", &
330  "A tool to create a MOLOPT or ADMM basis for a given set"// &
331  " of training structures", &
332  "Runs Tree Monte Carlo algorithm using additional input file(s)", &
333  "Runs (Tree) Monte Carlo trajectory file analysis", &
334  "Runs swarm based calculation"), &
335  enum_i_vals=(/do_atom, do_farming, do_test, do_cp2k, do_optimize_input, &
337  default_i_val=do_cp2k)
338  CALL section_add_keyword(section, keyword)
339  CALL keyword_release(keyword)
340 
341  CALL keyword_create(keyword, __location__, name="PROJECT_NAME", &
342  variants=(/"PROJECT"/), &
343  description="Name of the project (used to build the name of the "// &
344  "trajectory, and other files generated by the program)", &
345  usage="PROJECT_NAME {STRING}", &
346  default_c_val="PROJECT")
347  CALL section_add_keyword(section, keyword)
348  CALL keyword_release(keyword)
349 
350  CALL keyword_create(keyword, __location__, name="OUTPUT_FILE_NAME", &
351  description="Name of the output file. "// &
352  "Relevant only if automatically started (through farming for example). "// &
353  "If empty uses the project name as basis for it.", &
354  usage="OUTPUT_FILE_NAME {filename}", default_lc_val="")
355  CALL section_add_keyword(section, keyword)
356  CALL keyword_release(keyword)
357 
358  CALL keyword_create( &
359  keyword, __location__, name="RUN_TYPE", &
360  description="Type of run that you want to perform Geometry "// &
361  "optimization, md, montecarlo,...", &
362  usage="RUN_TYPE MD", &
363  default_i_val=energy_force_run, &
364  citations=(/ceriotti2014, schonherr2014/), &
365  enum_c_vals=s2a("NONE", "ENERGY", "ENERGY_FORCE", "MD", "GEO_OPT", &
366  "MC", "SPECTRA", "DEBUG", "BSSE", "LR", "PINT", "VIBRATIONAL_ANALYSIS", &
367  "BAND", "CELL_OPT", "WFN_OPT", "WAVEFUNCTION_OPTIMIZATION", &
368  "MOLECULAR_DYNAMICS", "GEOMETRY_OPTIMIZATION", "MONTECARLO", &
369  "ELECTRONIC_SPECTRA", "LINEAR_RESPONSE", "NORMAL_MODES", "RT_PROPAGATION", &
370  "EHRENFEST_DYN", "TAMC", "TMC", "DRIVER", "NEGF"), &
371  enum_i_vals=(/none_run, energy_run, energy_force_run, mol_dyn_run, &
377  enum_desc=s2a("Perform no tasks", "Computes energy", "Computes energy and forces", &
378  "Molecular Dynamics", "Geometry Optimization", "Monte Carlo", "Computes absorption Spectra", &
379  "Performs a Debug analysis", "Basis set superposition error", "Linear Response", &
380  "Path integral", "Vibrational analysis", "Band methods", &
381  "Cell optimization. Both cell vectors and atomic positions are optimised.", &
382  "Alias for ENERGY", "Alias for ENERGY", "Alias for MD", "Alias for GEO_OPT", &
383  "Alias for MC", "Alias for SPECTRA", "Alias for LR", "Alias for VIBRATIONAL_ANALYSIS", &
384  "Real Time propagation run (fixed ionic positions)", &
385  "Ehrenfest dynamics (using real time propagation of the wavefunction)", &
386  "Temperature Accelerated Monte Carlo (TAMC)", &
387  "Tree Monte Carlo (TMC), a pre-sampling MC algorithm", &
388  "i-PI driver mode", &
389  "Non-equilibrium Green's function method"))
390  CALL section_add_keyword(section, keyword)
391  CALL keyword_release(keyword)
392 
393  CALL keyword_create(keyword, __location__, name="WALLTIME", &
394  variants=(/"WALLTI"/), &
395  description="Maximum execution time for this run. Time in seconds or in HH:MM:SS.", &
396  usage="WALLTIME {real} or {HH:MM:SS}", default_lc_val="")
397  CALL section_add_keyword(section, keyword)
398  CALL keyword_release(keyword)
399 
400  CALL keyword_create(keyword, __location__, name="ECHO_INPUT", &
401  description="If the input should be echoed to the output with all the "// &
402  "defaults made explicit", &
403  usage="ECHO_INPUT NO", default_l_val=.false., lone_keyword_l_val=.true.)
404  CALL section_add_keyword(section, keyword)
405  CALL keyword_release(keyword)
406 
407  CALL keyword_create(keyword, __location__, name="ECHO_ALL_HOSTS", &
408  description="Echo a list of hostname and pid for all MPI processes.", &
409  usage="ECHO_ALL_HOSTS NO", default_l_val=.false., lone_keyword_l_val=.true.)
410  CALL section_add_keyword(section, keyword)
411  CALL keyword_release(keyword)
412 
413  CALL keyword_create(keyword, __location__, name="ENABLE_MPI_IO", &
414  description="Enable MPI parallelization for all supported I/O routines "// &
415  "Currently, only cube file writer/reader routines use MPI I/O. Disabling "// &
416  "this flag might speed up calculations dominated by I/O.", &
417  usage="ENABLE_MPI_IO FALSE", default_l_val=.true., lone_keyword_l_val=.true.)
418  CALL section_add_keyword(section, keyword)
419  CALL keyword_release(keyword)
420 
421  CALL keyword_create(keyword, __location__, name="TRACE", &
422  description="If a debug trace of the execution of the program should be written ", &
423  usage="TRACE", &
424  default_l_val=.false., lone_keyword_l_val=.true.)
425  CALL section_add_keyword(section, keyword)
426  CALL keyword_release(keyword)
427 
428  CALL keyword_create(keyword, __location__, name="TRACE_MASTER", &
429  description="For parallel TRACEd runs: only the master node writes output.", &
430  usage="TRACE_MASTER", &
431  default_l_val=.true., lone_keyword_l_val=.true.)
432  CALL section_add_keyword(section, keyword)
433  CALL keyword_release(keyword)
434 
435  CALL keyword_create( &
436  keyword, __location__, name="TRACE_MAX", &
437  description="Limit the total number a given subroutine is printed in the trace. Accounting is not influenced.", &
438  usage="TRACE_MAX 100", default_i_val=huge(0))
439  CALL section_add_keyword(section, keyword)
440  CALL keyword_release(keyword)
441 
442  CALL keyword_create( &
443  keyword, __location__, name="TRACE_ROUTINES", &
444  description="A list of routines to trace. If left empty all routines are traced. Accounting is not influenced.", &
445  usage="TRACE_ROUTINES {routine_name1} {routine_name2} ...", type_of_var=char_t, &
446  n_var=-1)
447  CALL section_add_keyword(section, keyword)
448  CALL keyword_release(keyword)
449 
450  CALL keyword_create( &
451  keyword, __location__, name="FLUSH_SHOULD_FLUSH", &
452  description="Flush output regularly, enabling this option might degrade performance significantly on certain machines.", &
453  usage="FLUSH_SHOULD_FLUSH", &
454  default_l_val=.true., lone_keyword_l_val=.true.)
455  CALL section_add_keyword(section, keyword)
456  CALL keyword_release(keyword)
457 
458  CALL keyword_create(keyword, __location__, name="CALLGRAPH", &
459  description="At the end of the run write a callgraph to file, "// &
460  "which contains detailed timing informations. "// &
461  "This callgraph can be viewed e.g. with the open-source program kcachegrind.", &
462  usage="CALLGRAPH {NONE|MASTER|ALL}", &
463  default_i_val=callgraph_none, lone_keyword_i_val=callgraph_master, &
464  enum_c_vals=s2a("NONE", "MASTER", "ALL"), &
465  enum_desc=s2a("No callgraph gets written", &
466  "Only the master process writes his callgraph", &
467  "All processes write their callgraph (into a separate files)."), &
469  CALL section_add_keyword(section, keyword)
470  CALL keyword_release(keyword)
471 
472  CALL keyword_create(keyword, __location__, name="CALLGRAPH_FILE_NAME", &
473  description="Name of the callgraph file, which is written at the end of the run. "// &
474  "If not specified the project name will be used as filename.", &
475  usage="CALLGRAPH_FILE_NAME {filename}", default_lc_val="")
476  CALL section_add_keyword(section, keyword)
477  CALL keyword_release(keyword)
478 
479  CALL keyword_create(keyword, __location__, name="SEED", &
480  description="Initial seed for the global (pseudo)random number generator "// &
481  "to create a stream of normally Gaussian distributed random numbers. "// &
482  "Exactly 1 or 6 positive integer values are expected. A single value is "// &
483  "replicated to fill up the full seed array with 6 numbers.", &
484  n_var=-1, &
485  type_of_var=integer_t, &
486  usage="SEED {INTEGER} .. {INTEGER}", &
487  default_i_vals=(/2000/))
488  CALL section_add_keyword(section, keyword)
489  CALL keyword_release(keyword)
490 
491  CALL keyword_create(keyword, __location__, name="SAVE_MEM", &
492  description="Some sections of the input structure are deallocated when not needed,"// &
493  " and reallocated only when used. This reduces the required maximum memory ", &
494  usage="SAVE_MEM", &
495  default_l_val=.false., lone_keyword_l_val=.true.)
496  CALL section_add_keyword(section, keyword)
497  CALL keyword_release(keyword)
498 
499  CALL cp_print_key_section_create(print_key, __location__, "TIMINGS", description= &
500  "Controls the printing of the timing report at the end of CP2K execution", &
501  print_level=silent_print_level, filename="__STD_OUT__")
502 
503  CALL keyword_create(keyword, __location__, name="THRESHOLD", &
504  description="Specify % of CPUTIME above which the contribution will be inserted in the"// &
505  " final timing report (e.g. 0.02 = 2%)", &
506  usage="THRESHOLD {REAL}", &
507  default_r_val=0.02_dp)
508  CALL section_add_keyword(print_key, keyword)
509  CALL keyword_release(keyword)
510 
511  CALL keyword_create(keyword, __location__, name="SORT_BY_SELF_TIME", &
512  description="Sort the final timing report by the average self (exclusive) time instead of the "// &
513  "total (inclusive) time of a routine", &
514  usage="SORT_BY_SELF_TIME on", &
515  default_l_val=.false., lone_keyword_l_val=.true.)
516  CALL section_add_keyword(print_key, keyword)
517  CALL keyword_release(keyword)
518 
519  CALL keyword_create(keyword, __location__, name="REPORT_MAXLOC", &
520  description="Report the rank with the slowest maximum self timing."// &
521  " Can be used to debug hard- or software."// &
522  " Also enables ECHO_ALL_HOSTS to link rank to hostname.", &
523  usage="REPORT_MAXLOC on", &
524  default_l_val=.false., lone_keyword_l_val=.true.)
525  CALL section_add_keyword(print_key, keyword)
526  CALL keyword_release(keyword)
527 
528  CALL keyword_create(keyword, __location__, name="TIME_MPI", &
529  description="Include message_passing calls in the timing report (useful with CALLGRAPH).", &
530  usage="TIME_MPI .FALSE.", &
531  default_l_val=.true., lone_keyword_l_val=.true.)
532  CALL section_add_keyword(print_key, keyword)
533  CALL keyword_release(keyword)
534 
535  CALL keyword_create(keyword, __location__, name="TIMINGS_LEVEL", &
536  description="Specify the level of timings report. "// &
537  "Possible values are: 0 (report only CP2K root timer), 1 (all timers).", &
538  usage="TIMINGS_LEVEL 1", &
539  default_i_val=default_timings_level, lone_keyword_i_val=default_timings_level)
540  CALL section_add_keyword(print_key, keyword)
541  CALL keyword_release(keyword)
542 
543  CALL section_add_subsection(section, print_key)
544  CALL section_release(print_key)
545 
546  CALL cp_print_key_section_create(print_key, __location__, "REFERENCES", description= &
547  "Controls the printing of the references relevant to the calculations performed", &
548  print_level=silent_print_level, filename="__STD_OUT__")
549  CALL section_add_subsection(section, print_key)
550  CALL section_release(print_key)
551 
552  CALL cp_print_key_section_create(print_key, __location__, "PROGRAM_RUN_INFO", &
553  description="controls the printing of initialization controlled by the global section", &
554  print_level=silent_print_level, add_last=add_last_numeric, filename="__STD_OUT__")
555  CALL section_add_subsection(section, print_key)
556  CALL section_release(print_key)
557 
558  CALL cp_print_key_section_create(print_key, __location__, "PRINT", description= &
559  "controls the printing of physical and mathematical constants", &
560  print_level=medium_print_level, filename="__STD_OUT__")
561 
562  CALL keyword_create(keyword, __location__, name="BASIC_DATA_TYPES", &
563  description="Controls the printing of the basic data types.", &
564  default_l_val=.false., lone_keyword_l_val=.true.)
565  CALL section_add_keyword(print_key, keyword)
566  CALL keyword_release(keyword)
567  CALL keyword_create(keyword, __location__, name="physcon", &
568  description="if the printkey is active prints the physical constants", &
569  default_l_val=.true., lone_keyword_l_val=.true.)
570  CALL section_add_keyword(print_key, keyword)
571  CALL keyword_release(keyword)
572  CALL keyword_create(keyword, __location__, name="SPHERICAL_HARMONICS", &
573  description="if the printkey is active prints the spherical harmonics", &
574  default_i_val=-1, type_of_var=integer_t)
575  CALL section_add_keyword(print_key, keyword)
576  CALL keyword_release(keyword)
577  CALL keyword_create(keyword, __location__, name="RNG_MATRICES", &
578  description="Prints the transformation matrices used by the random number generator", &
579  default_l_val=.false., &
580  lone_keyword_l_val=.true.)
581  CALL section_add_keyword(print_key, keyword)
582  CALL keyword_release(keyword)
583  CALL keyword_create(keyword, __location__, name="RNG_CHECK", &
584  description="Performs a check of the global (pseudo)random "// &
585  "number generator (RNG) and prints the result", &
586  default_l_val=.false., &
587  lone_keyword_l_val=.true.)
588  CALL section_add_keyword(print_key, keyword)
589  CALL keyword_release(keyword)
590  CALL keyword_create(keyword, __location__, name="GLOBAL_GAUSSIAN_RNG", &
591  description="Prints the initial status of the global Gaussian "// &
592  "(pseudo)random number stream which is mostly used for "// &
593  "the velocity initialization", &
594  default_l_val=.false., &
595  lone_keyword_l_val=.true.)
596  CALL section_add_keyword(print_key, keyword)
597  CALL keyword_release(keyword)
598 
599  CALL section_add_subsection(section, print_key)
600  CALL section_release(print_key)
601  NULLIFY (sub_section)
602  ! FM section
603  CALL create_fm_section(sub_section)
604  CALL section_add_subsection(section, sub_section)
605  CALL section_release(sub_section)
606  ! DBCSR options
607  CALL create_dbcsr_section(sub_section)
608  CALL section_add_subsection(section, sub_section)
609  CALL section_release(sub_section)
610  ! FM diagonalization redistribution rules
611  CALL create_fm_diag_rules_section(sub_section)
612  CALL section_add_subsection(section, sub_section)
613  CALL section_release(sub_section)
614  ! Grid library
615  CALL create_grid_section(sub_section)
616  CALL section_add_subsection(section, sub_section)
617  CALL section_release(sub_section)
618 
619  END SUBROUTINE create_global_section
620 
621 ! **************************************************************************************************
622 !> \brief Creates the dbcsr section for configuring FM
623 !> \param section ...
624 !> \date 2011-04-05
625 !> \author Florian Schiffmann
626 ! **************************************************************************************************
627  SUBROUTINE create_fm_section(section)
628  TYPE(section_type), POINTER :: section
629 
630  INTEGER :: default_matmul
631  TYPE(keyword_type), POINTER :: keyword
632 
633  cpassert(.NOT. ASSOCIATED(section))
634  CALL section_create(section, __location__, name="FM", &
635  description="Configuration options for the full matrices.", &
636  n_keywords=1, n_subsections=0, repeats=.false.)
637 
638  NULLIFY (keyword)
639 
640  CALL keyword_create(keyword, __location__, name="NROW_BLOCKS", &
641  description="Defines the number of rows per scalapack block in "// &
642  "the creation of block cyclic dense matrices ", &
643  default_i_val=64)
644  CALL section_add_keyword(section, keyword)
645  CALL keyword_release(keyword)
646 
647  CALL keyword_create(keyword, __location__, name="NCOL_BLOCKS", &
648  description="Defines the number of columns per scalapack block in "// &
649  "the creation of vlock cyclic dense matrices ", &
650  default_i_val=64)
651  CALL section_add_keyword(section, keyword)
652  CALL keyword_release(keyword)
653 
654  CALL keyword_create(keyword, __location__, name="FORCE_BLOCK_SIZE", &
655  description="Ensure for small matrices that the layout is compatible "// &
656  "with bigger ones, i.e. no subdivision is performed (can break LAPACK!!!).", &
657  usage="FORCE_BLOCK_SIZE", &
658  default_l_val=.false., lone_keyword_l_val=.true.)
659  CALL section_add_keyword(section, keyword)
660  CALL keyword_release(keyword)
661 
662 #if defined(__COSMA)
663  default_matmul = do_cosma
664 #else
665  default_matmul = do_scalapack
666 #endif
667 
668  CALL keyword_create(keyword, __location__, name="TYPE_OF_MATRIX_MULTIPLICATION", &
669  description="Allows to switch between scalapack pxgemm and COSMA pxgemm. "// &
670  "COSMA reduces the communication costs but increases the memory demands. "// &
671  "The performance of Scalapack's pxgemm on GPU's depends "// &
672  "crucially on the BLOCK_SIZES. Make sure optimized kernels are available.", &
673  default_i_val=default_matmul, &
674  enum_i_vals=(/do_scalapack, do_scalapack, do_cosma/), &
675  enum_c_vals=s2a("SCALAPACK", "PDGEMM", "COSMA"), &
676  enum_desc=s2a("Standard ScaLAPACK pdgemm", &
677  "Alias for ScaLAPACK", &
678  "COSMA is employed. See <https://github.com/eth-cscs/COSMA>."))
679  CALL section_add_keyword(section, keyword)
680  CALL keyword_release(keyword)
681 
682  !
683  END SUBROUTINE create_fm_section
684 ! **************************************************************************************************
685 !> \brief Creates the input section used to define the heuristic rules which determine if
686 !> a FM matrix should be redistributed before diagonalizing it.
687 !> \param section the input section to create
688 !> \author Nico Holmberg [01.2018]
689 ! **************************************************************************************************
690  SUBROUTINE create_fm_diag_rules_section(section)
691  TYPE(section_type), POINTER :: section
692 
693  TYPE(keyword_type), POINTER :: keyword
694 
695  cpassert(.NOT. ASSOCIATED(section))
696  CALL section_create(section, __location__, name="FM_DIAG_SETTINGS", &
697  description="This section defines a set of heuristic rules which are "// &
698  "used to calculate the optimal number of CPUs, M, needed to diagonalize a "// &
699  "full matrix distributed on N processors (FM type). If M &lt N, the matrix "// &
700  "is redistributed onto M processors before it is diagonalized. "// &
701  "The optimal value is calculate according to M = ((K+a*x-1)/(a*x))*a, "// &
702  "where K is the size of the matrix, and {a, x} are integers defined below. "// &
703  "The default values have been selected based on timings on a Cray XE6. "// &
704  "Supports diagonalization libraries SL and ELPA (see keyword ELPA_FORCE_REDISTRIBUTE).", &
705  n_keywords=3, n_subsections=0, repeats=.false.)
706 
707  NULLIFY (keyword)
708 
709  CALL keyword_create(keyword, __location__, name="PARAMETER_A", &
710  description="Parameter used for defining the rule which determines the optimal "// &
711  "number of CPUs needed to diagonalize a full distributed matrix. The optimal "// &
712  "number of CPUs will be an integer multiple of this variable.", &
713  usage="PARAMETER_A 4", type_of_var=integer_t, &
714  default_i_val=4)
715  CALL section_add_keyword(section, keyword)
716  CALL keyword_release(keyword)
717 
718  CALL keyword_create(keyword, __location__, name="PARAMETER_X", &
719  description="Parameter used for defining the rule which determines the optimal "// &
720  "number of CPUs needed to diagonalize a full distributed matrix. The optimal "// &
721  "number of CPUs will be roughly proportional to this value.", &
722  usage="PARAMETER_X 60", type_of_var=integer_t, &
723  default_i_val=60)
724  CALL section_add_keyword(section, keyword)
725  CALL keyword_release(keyword)
726 
727  CALL keyword_create(keyword, __location__, name="PRINT_FM_REDISTRIBUTE", &
728  description="Controls printing of information related to this section. For each "// &
729  "diagonalized matrix, prints the size of the matrix, the optimal number of CPUs, "// &
730  "as well as notifies if the matrix was redistributed. Useful for testing.", &
731  usage="PRINT_FM_REDISTRIBUTE", type_of_var=logical_t, &
732  default_l_val=.false., lone_keyword_l_val=.true.)
733  CALL section_add_keyword(section, keyword)
734  CALL keyword_release(keyword)
735 
736  CALL keyword_create(keyword, __location__, name="ELPA_FORCE_REDISTRIBUTE", &
737  description="Controls how to perform redistribution when ELPA is used for diagonalization. "// &
738  "By default, redistribution is always performed using the defined rules. "// &
739  "By turning off this keyword, matrices are redistributed only to prevent crashes in the ELPA "// &
740  "library which happens when the original matrix is distributed over too many processors. ", &
741  usage="ELPA_FORCE_REDISTRIBUTE", type_of_var=logical_t, &
742  default_l_val=.true., lone_keyword_l_val=.true.)
743  CALL section_add_keyword(section, keyword)
744  CALL keyword_release(keyword)
745 
746  END SUBROUTINE create_fm_diag_rules_section
747 
748 ! **************************************************************************************************
749 !> \brief Creates the section for configuring the grid library
750 !> \param section ...
751 !> \author Ole Schuett
752 ! **************************************************************************************************
753  SUBROUTINE create_grid_section(section)
754  TYPE(section_type), POINTER :: section
755 
756  TYPE(keyword_type), POINTER :: keyword
757 
758  cpassert(.NOT. ASSOCIATED(section))
759  CALL section_create(section, __location__, name="GRID", &
760  description="Configuration options for the grid library, "// &
761  "which performs e.g. the collocate and integrate of the GPW method.", &
762  n_keywords=1, n_subsections=0, repeats=.false.)
763 
764  NULLIFY (keyword)
765  CALL keyword_create(keyword, __location__, name="BACKEND", &
766  description="Selects the backed used by the grid library.", &
767  default_i_val=grid_backend_auto, &
770  enum_c_vals=s2a("AUTO", "REFERENCE", "CPU", "DGEMM", "GPU", "HIP"), &
771  enum_desc=s2a("Let the grid library pick the backend automatically", &
772  "Reference backend implementation", &
773  "Optimized CPU backend", &
774  "Alternative CPU backend based on DGEMM", &
775  "GPU backend optimized for CUDA that also supports HIP", &
776  "HIP backend optimized for ROCm"))
777  CALL section_add_keyword(section, keyword)
778  CALL keyword_release(keyword)
779 
780  CALL keyword_create(keyword, __location__, name="VALIDATE", &
781  description="When enabled the reference backend is run in shadow mode "// &
782  "and its results are compared with those from the selected backend. "// &
783  "If the two results differ by too much then the calculation is aborted.", &
784  default_l_val=.false., lone_keyword_l_val=.true.)
785  CALL section_add_keyword(section, keyword)
786  CALL keyword_release(keyword)
787 
788  CALL keyword_create(keyword, __location__, name="APPLY_CUTOFF", &
789  description="When enabled the cpu backend "// &
790  "apply a spherical cutoff on the top of the cube. "// &
791  "There is a performance penalty using it in "// &
792  "combination with the cpu backend but it is on by "// &
793  "default for the regtests", default_l_val=.true., &
794  lone_keyword_l_val=.true.)
795  CALL section_add_keyword(section, keyword)
796  CALL keyword_release(keyword)
797 
798  END SUBROUTINE create_grid_section
799 
800 END MODULE input_cp2k_global
collects all references to literature in CP2K as new algorithms / method are included from literature...
Definition: bibliography.F:28
integer, save, public schonherr2014
Definition: bibliography.F:43
integer, save, public frigo2005
Definition: bibliography.F:43
integer, save, public ceriotti2014
Definition: bibliography.F:43
methods related to the blacs parallel environment
Definition: cp_blacs_env.F:15
integer, parameter, public blacs_grid_row
Definition: cp_blacs_env.F:32
integer, parameter, public blacs_grid_col
Definition: cp_blacs_env.F:32
integer, parameter, public blacs_grid_square
Definition: cp_blacs_env.F:32
used for collecting some of the diagonalization schemes available for cp_fm_type. cp_fm_power also mo...
Definition: cp_fm_diag.F:17
integer, parameter, public fm_diag_type_cusolver
Definition: cp_fm_diag.F:87
integer, parameter, public fm_diag_type_dlaf
Definition: cp_fm_diag.F:87
integer, parameter, public fm_diag_type_scalapack
Definition: cp_fm_diag.F:87
integer, parameter, public fm_diag_type_default
Definition: cp_fm_diag.F:98
integer, parameter, public fm_diag_type_elpa
Definition: cp_fm_diag.F:87
Wrapper for ELPA.
Definition: cp_fm_elpa.F:12
character(len=14), dimension(1), parameter, public elpa_kernel_names
Definition: cp_fm_elpa.F:150
character(len=44), dimension(1), parameter, public elpa_kernel_descriptions
Definition: cp_fm_elpa.F:151
integer, dimension(1), parameter, public elpa_kernel_ids
Definition: cp_fm_elpa.F:149
routines to handle the output, The idea is to remove the decision of wheter to output and what to out...
integer, parameter, public debug_print_level
integer, parameter, public low_print_level
integer, parameter, public medium_print_level
integer, parameter, public high_print_level
integer, parameter, public add_last_numeric
integer, parameter, public silent_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
Fortran API for the grid package, which is written in C.
Definition: grid_api.F:12
integer, parameter, public grid_backend_auto
Definition: grid_api.F:66
integer, parameter, public grid_backend_gpu
Definition: grid_api.F:70
integer, parameter, public grid_backend_hip
Definition: grid_api.F:71
integer, parameter, public grid_backend_dgemm
Definition: grid_api.F:69
integer, parameter, public grid_backend_cpu
Definition: grid_api.F:68
integer, parameter, public grid_backend_ref
Definition: grid_api.F:67
collects all constants needed in input so that they can be used without circular dependencies
integer, parameter, public fftw_plan_patient
integer, parameter, public driver_run
integer, parameter, public energy_run
integer, parameter, public fftw_plan_exhaustive
integer, parameter, public do_fft_sg
integer, parameter, public callgraph_all
integer, parameter, public do_farming
integer, parameter, public electronic_spectra_run
integer, parameter, public do_cosma
integer, parameter, public do_cp2k
integer, parameter, public do_scalapack
integer, parameter, public linear_response_run
integer, parameter, public do_tamc
integer, parameter, public ehrenfest
integer, parameter, public do_opt_basis
integer, parameter, public debug_run
integer, parameter, public do_test
integer, parameter, public fftw_plan_estimate
integer, parameter, public callgraph_master
integer, parameter, public do_dgemm_blas
integer, parameter, public do_band
integer, parameter, public do_tree_mc
integer, parameter, public bsse_run
integer, parameter, public energy_force_run
integer, parameter, public callgraph_none
integer, parameter, public do_optimize_input
integer, parameter, public do_atom
integer, parameter, public mon_car_run
integer, parameter, public do_tree_mc_ana
integer, parameter, public mol_dyn_run
integer, parameter, public do_fft_fftw3
integer, parameter, public cell_opt_run
integer, parameter, public pint_run
integer, parameter, public gaussian
integer, parameter, public vib_anal
integer, parameter, public none_run
integer, parameter, public negf_run
integer, parameter, public fftw_plan_measure
integer, parameter, public tree_mc_run
integer, parameter, public real_time_propagation
integer, parameter, public geo_opt_run
integer, parameter, public do_dgemm_spla
integer, parameter, public do_swarm
builds the global input section for cp2k
subroutine, public create_global_section(section)
section to hold global settings for the whole program
represents keywords in an input
subroutine, public keyword_release(keyword)
releases the given keyword (see doc/ReferenceCounting.html)
subroutine, public keyword_create(keyword, location, name, description, usage, type_of_var, n_var, repeats, variants, default_val, default_l_val, default_r_val, default_lc_val, default_c_val, default_i_val, default_l_vals, default_r_vals, default_c_vals, default_i_vals, lone_keyword_val, lone_keyword_l_val, lone_keyword_r_val, lone_keyword_c_val, lone_keyword_i_val, lone_keyword_l_vals, lone_keyword_r_vals, lone_keyword_c_vals, lone_keyword_i_vals, enum_c_vals, enum_i_vals, enum, enum_strict, enum_desc, unit_str, citations, deprecation_notice, removed)
creates a keyword object
objects that represent the structure of input sections and the data contained in an input section
subroutine, public section_create(section, location, name, description, n_keywords, n_subsections, repeats, citations)
creates a list of keywords
subroutine, public section_add_keyword(section, keyword)
adds a keyword to the given section
subroutine, public section_add_subsection(section, subsection)
adds a subsection to the given section
recursive subroutine, public section_release(section)
releases the given keyword list (see doc/ReferenceCounting.html)
a wrapper for basic fortran types.
integer, parameter, public logical_t
integer, parameter, public char_t
integer, parameter, public integer_t
Defines the basic variable types.
Definition: kinds.F:23
integer, parameter, public dp
Definition: kinds.F:34
Utilities for string manipulations.
Timing routines for accounting.
Definition: timings.F:17
integer, parameter, public default_timings_level
Definition: timings.F:67