(git:374b731)
Loading...
Searching...
No Matches
cp_dbcsr_cp2k_link.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 Routines that link DBCSR and CP2K concepts together
10!> \author Ole Schuett
11!> \par History
12!> 01.2014 created
13! **************************************************************************************************
15 USE ao_util, ONLY: exp_radius
20 USE bibliography, ONLY: borstnik2014,&
23 cite_reference
26 USE dbcsr_api, ONLY: &
27 dbcsr_finalize, dbcsr_get_block_p, dbcsr_get_default_config, dbcsr_get_matrix_type, &
28 dbcsr_has_symmetry, dbcsr_reserve_blocks, dbcsr_set, dbcsr_set_config, dbcsr_type, &
29 dbcsr_type_no_symmetry
41 USE kinds, ONLY: dp,&
42 real_8
43 USE orbital_pointers, ONLY: nso
46 USE qs_ks_types, ONLY: get_ks_env,&
55 USE string_utilities, ONLY: s2a
56#include "./base/base_uses.f90"
57
58 IMPLICIT NONE
59
60 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'cp_dbcsr_cp2k_link'
61
62 PUBLIC :: create_dbcsr_section
63 PUBLIC :: cp_dbcsr_config
66
67 PRIVATE
68
69 ! Possible drivers to use for matrix multiplications
70 INTEGER, PARAMETER :: mm_driver_auto = 0
71 INTEGER, PARAMETER :: mm_driver_matmul = 1
72 INTEGER, PARAMETER :: mm_driver_blas = 2
73 INTEGER, PARAMETER :: mm_driver_smm = 3
74 INTEGER, PARAMETER :: mm_driver_xsmm = 4
75
76 CHARACTER(len=*), PARAMETER :: mm_name_auto = "AUTO", &
77 mm_name_blas = "BLAS", &
78 mm_name_matmul = "MATMUL", &
79 mm_name_smm = "SMM", &
80 mm_name_xsmm = "XSMM"
81CONTAINS
82
83! **************************************************************************************************
84!> \brief Creates the dbcsr section for configuring DBCSR
85!> \param section ...
86!> \date 2011-04-05
87!> \author Urban Borstnik
88! **************************************************************************************************
89 SUBROUTINE create_dbcsr_section(section)
90 TYPE(section_type), POINTER :: section
91
92 INTEGER :: idefault
93 LOGICAL :: ldefault
94 REAL(kind=dp) :: rdefault
95 TYPE(keyword_type), POINTER :: keyword
96 TYPE(section_type), POINTER :: subsection
97
98 cpassert(.NOT. ASSOCIATED(section))
99 CALL section_create(section, __location__, name="DBCSR", &
100 description="Configuration options for the DBCSR library.", &
101 n_keywords=1, n_subsections=0, repeats=.false., &
102 citations=(/borstnik2014, schuett2016/))
103
104 NULLIFY (keyword)
105
106 CALL keyword_create(keyword, __location__, name="mm_stack_size", &
107 description="Size of multiplication parameter stack." &
108 //" A negative value leaves the decision up to DBCSR.", &
109 usage="mm_stack_size 1000", &
110 default_i_val=-1)
111 CALL section_add_keyword(section, keyword)
112 CALL keyword_release(keyword)
113
114 CALL keyword_create(keyword, __location__, name="mm_driver", &
115 description="Select which backend to use preferably "// &
116 "for matrix block multiplications on the host.", &
117 usage="mm_driver blas", &
118 default_i_val=mm_driver_auto, &
119 enum_c_vals=s2a("AUTO", "BLAS", "MATMUL", "SMM", "XSMM"), &
120 enum_i_vals=(/mm_driver_auto, mm_driver_blas, mm_driver_matmul, mm_driver_smm, &
121 mm_driver_xsmm/), &
122 enum_desc=s2a("Choose automatically the best available driver", &
123 "BLAS (requires the BLAS library at link time)", &
124 "Fortran MATMUL", &
125 "Library optimised for Small Matrix Multiplies "// &
126 "(requires the SMM library at link time)", &
127 "LIBXSMM"), &
128 citations=(/heinecke2016/))
129 CALL section_add_keyword(section, keyword)
130 CALL keyword_release(keyword)
131
132 CALL dbcsr_get_default_config(avg_elements_images=idefault)
133 CALL keyword_create(keyword, __location__, name="avg_elements_images", &
134 description="Average number of elements (dense limit)" &
135 //" for each image, which also corresponds to" &
136 //" the average number of elements exchanged between MPI processes" &
137 //" during the operations." &
138 //" A negative or zero value means unlimited.", &
139 usage="avg_elements_images 10000", &
140 default_i_val=idefault)
141 CALL section_add_keyword(section, keyword)
142 CALL keyword_release(keyword)
143
144 CALL dbcsr_get_default_config(num_mult_images=idefault)
145 CALL keyword_create(keyword, __location__, name="num_mult_images", &
146 description="Multiplicative factor for number of virtual images.", &
147 usage="num_mult_images 2", &
148 default_i_val=idefault)
149 CALL section_add_keyword(section, keyword)
150 CALL keyword_release(keyword)
151
152 CALL dbcsr_get_default_config(use_mpi_allocator=ldefault)
153 CALL keyword_create(keyword, __location__, name="use_mpi_allocator", &
154 description="Use MPI allocator" &
155 //" to allocate buffers used in MPI communications.", &
156 usage="use_mpi_allocator T", &
157 default_l_val=ldefault)
158 CALL section_add_keyword(section, keyword)
159 CALL keyword_release(keyword)
160
161 CALL dbcsr_get_default_config(use_mpi_rma=ldefault)
162 CALL keyword_create(keyword, __location__, name="use_mpi_rma", &
163 description="Use RMA for MPI communications" &
164 //" for each image, which also corresponds to" &
165 //" the number of elements exchanged between MPI processes" &
166 //" during the operations.", &
167 usage="use_mpi_rma F", &
168 default_l_val=ldefault)
169 CALL section_add_keyword(section, keyword)
170 CALL keyword_release(keyword)
171
172 CALL dbcsr_get_default_config(num_layers_3d=idefault)
173 CALL keyword_create(keyword, __location__, name="num_layers_3D", &
174 description="Number of layers for the 3D multplication algorithm.", &
175 usage="num_layers_3D 1", &
176 default_i_val=idefault)
177 CALL section_add_keyword(section, keyword)
178 CALL keyword_release(keyword)
179
180 CALL dbcsr_get_default_config(nstacks=idefault)
181 CALL keyword_create(keyword, __location__, name="n_size_mnk_stacks", &
182 description="Number of stacks to use for distinct atomic sizes" &
183 //" (e.g., 2 for a system of mostly waters). ", &
184 usage="n_size_mnk_stacks 2", &
185 default_i_val=idefault)
186 CALL section_add_keyword(section, keyword)
187 CALL keyword_release(keyword)
188
189 CALL dbcsr_get_default_config(use_comm_thread=ldefault)
190 CALL keyword_create(keyword, __location__, name="use_comm_thread", &
191 description="During multiplication, use a thread to periodically poll" &
192 //" MPI to progress outstanding message completions. This is" &
193 //" beneficial on systems without a DMA-capable network adapter" &
194 //" e.g. Cray XE6.", &
195 usage="use_comm_thread T", &
196 default_l_val=ldefault)
197 CALL section_add_keyword(section, keyword)
198 CALL keyword_release(keyword)
199
200 CALL keyword_create(keyword, __location__, name="MAX_ELEMENTS_PER_BLOCK", &
201 description="Default block size for turning dense matrices in blocked ones", &
202 usage="MAX_ELEMENTS_PER_BLOCK 32", &
203 default_i_val=max_elements_per_block)
204 CALL section_add_keyword(section, keyword)
205 CALL keyword_release(keyword)
206
207 CALL keyword_create(keyword, __location__, name="comm_thread_load", &
208 description="If a communications thread is used, specify how much " &
209 //"multiplication workload (%) the thread should perform in " &
210 //"addition to communication tasks. " &
211 //"A negative value leaves the decision up to DBCSR.", &
212 usage="comm_thread_load 50", &
213 default_i_val=-1)
214 CALL section_add_keyword(section, keyword)
215 CALL keyword_release(keyword)
216
217 CALL dbcsr_get_default_config(multrec_limit=idefault)
218 CALL keyword_create(keyword, __location__, name="multrec_limit", &
219 description="Recursion limit of cache oblivious multrec algorithm.", &
220 default_i_val=idefault)
221 CALL section_add_keyword(section, keyword)
222 CALL keyword_release(keyword)
223
224 CALL dbcsr_get_default_config(use_mempools_cpu=ldefault)
225 CALL keyword_create(keyword, __location__, name="use_mempools_cpu", &
226 description="Enable memory pools on the CPU.", &
227 default_l_val=ldefault)
228 CALL section_add_keyword(section, keyword)
229 CALL keyword_release(keyword)
230
231 NULLIFY (subsection)
232 CALL section_create(subsection, __location__, name="TENSOR", &
233 description="Configuration options for Tensors.", &
234 n_keywords=1, n_subsections=0, repeats=.false.)
235
236 CALL dbcsr_get_default_config(tas_split_factor=rdefault)
237 CALL keyword_create(keyword, __location__, name="TAS_SPLIT_FACTOR", &
238 description="Parameter for hybrid DBCSR-TAS matrix-matrix multiplication algorithm: "// &
239 "a TAS matrix is split into s submatrices with s = N_max/(N_min*f) with f "// &
240 "given by this parameter and N_max/N_min the max/min occupancies of the matrices "// &
241 "involved in a multiplication. A large value makes the multiplication Cannon-based "// &
242 "(s=1) and a small value (> 0) makes the multiplication based on TAS algorithm "// &
243 "(s=number of MPI ranks)", &
244 default_r_val=rdefault)
245 CALL section_add_keyword(subsection, keyword)
246 CALL keyword_release(keyword)
247
248 CALL section_add_subsection(section, subsection)
249 CALL section_release(subsection)
250
251 !---------------------------------------------------------------------------
252 NULLIFY (subsection)
253 CALL section_create(subsection, __location__, name="ACC", &
254 description="Configuration options for the ACC-Driver.", &
255 n_keywords=1, n_subsections=0, repeats=.false.)
256
257 CALL dbcsr_get_default_config(accdrv_thread_buffers=idefault)
258 CALL keyword_create(keyword, __location__, name="thread_buffers", &
259 description="Number of transfer-buffers associated with each thread and corresponding stream.", &
260 default_i_val=idefault)
261 CALL section_add_keyword(subsection, keyword)
262 CALL keyword_release(keyword)
263
264 CALL dbcsr_get_default_config(accdrv_avoid_after_busy=ldefault)
265 CALL keyword_create(keyword, __location__, name="avoid_after_busy", &
266 description="If enabled, stacks are not processed by the acc-driver " &
267 //"after it has signaled congestion during a round of flushing. " &
268 //"For the next round of flusing the driver is used again.", &
269 default_l_val=ldefault)
270 CALL section_add_keyword(subsection, keyword)
271 CALL keyword_release(keyword)
272
273 CALL dbcsr_get_default_config(accdrv_min_flop_process=idefault)
274 CALL keyword_create(keyword, __location__, name="min_flop_process", &
275 description="Only process stacks with more than the given number of " &
276 //"floating-point operations per stack-entry (2*m*n*k).", &
277 default_i_val=idefault)
278 CALL section_add_keyword(subsection, keyword)
279 CALL keyword_release(keyword)
280
281 CALL dbcsr_get_default_config(accdrv_stack_sort=ldefault)
282 CALL keyword_create(keyword, __location__, name="stack_sort", &
283 description="Sort multiplication stacks according to C-access.", &
284 default_l_val=ldefault)
285 CALL section_add_keyword(subsection, keyword)
286 CALL keyword_release(keyword)
287
288 CALL dbcsr_get_default_config(accdrv_min_flop_sort=idefault)
289 CALL keyword_create(keyword, __location__, name="min_flop_sort", &
290 description="Only sort stacks with more than the given number of " &
291 //"floating-point operations per stack-entry (2*m*n*k). " &
292 //"Alternatively, the stacks are roughly ordered through a " &
293 //"binning-scheme by Peter Messmer. (Depends on ACC%STACK_SORT)", &
294 default_i_val=idefault)
295 CALL section_add_keyword(subsection, keyword)
296 CALL keyword_release(keyword)
297
298 CALL dbcsr_get_default_config(accdrv_do_inhomogenous=ldefault)
299 CALL keyword_create(keyword, __location__, name="process_inhomogenous", &
300 description="If enabled, inhomogenous stacks are also processed by the acc driver.", &
301 default_l_val=ldefault)
302 CALL section_add_keyword(subsection, keyword)
303 CALL keyword_release(keyword)
304
305 CALL dbcsr_get_default_config(accdrv_binning_nbins=idefault)
306 CALL keyword_create(keyword, __location__, name="binning_nbins", &
307 description="Number of bins used when ordering " &
308 //"the stacks with the binning-scheme.", &
309 default_i_val=idefault)
310 CALL section_add_keyword(subsection, keyword)
311 CALL keyword_release(keyword)
312
313 CALL dbcsr_get_default_config(accdrv_binning_binsize=idefault)
314 CALL keyword_create(keyword, __location__, name="binning_binsize", &
315 description="Size of bins used when ordering " &
316 //"the stacks with the binning-scheme.", &
317 default_i_val=idefault)
318 CALL section_add_keyword(subsection, keyword)
319 CALL keyword_release(keyword)
320
321 CALL section_add_subsection(section, subsection)
322 CALL section_release(subsection)
323
324 END SUBROUTINE create_dbcsr_section
325
326! **************************************************************************************************
327!> \brief Configures options for DBCSR
328!> \param root_section ...
329! **************************************************************************************************
330 SUBROUTINE cp_dbcsr_config(root_section)
331 TYPE(section_vals_type), POINTER :: root_section
332
333 CHARACTER(len=*), PARAMETER :: routinen = 'cp_dbcsr_config'
334
335 INTEGER :: handle, ival
336 LOGICAL :: lval
337 REAL(kind=real_8) :: rval
338 TYPE(section_vals_type), POINTER :: dbcsr_section
339
340 CALL timeset(routinen, handle)
341
342 CALL cite_reference(borstnik2014)
343 CALL cite_reference(schuett2016)
344
345 dbcsr_section => section_vals_get_subs_vals(root_section, "GLOBAL%DBCSR")
346
347 CALL section_vals_val_get(dbcsr_section, "mm_stack_size", i_val=ival)
348 CALL dbcsr_set_config(mm_stack_size=ival)
349
350 CALL section_vals_val_get(dbcsr_section, "MAX_ELEMENTS_PER_BLOCK", i_val=max_elements_per_block)
351
352 CALL section_vals_val_get(dbcsr_section, "avg_elements_images", i_val=ival)
353 CALL dbcsr_set_config(avg_elements_images=ival)
354
355 CALL section_vals_val_get(dbcsr_section, "num_mult_images", i_val=ival)
356 CALL dbcsr_set_config(num_mult_images=ival)
357
358 CALL section_vals_val_get(dbcsr_section, "n_size_mnk_stacks", i_val=ival)
359 CALL dbcsr_set_config(nstacks=ival)
360
361 CALL section_vals_val_get(dbcsr_section, "use_mpi_allocator", l_val=lval)
362 CALL dbcsr_set_config(use_mpi_allocator=lval)
363
364 CALL section_vals_val_get(dbcsr_section, "use_mpi_rma", l_val=lval)
365 CALL dbcsr_set_config(use_mpi_rma=lval)
366
367 CALL section_vals_val_get(dbcsr_section, "num_layers_3D", i_val=ival)
368 CALL dbcsr_set_config(num_layers_3d=ival)
369
370 CALL section_vals_val_get(dbcsr_section, "use_comm_thread", l_val=lval)
371 CALL dbcsr_set_config(use_comm_thread=lval)
372
373 CALL section_vals_val_get(dbcsr_section, "comm_thread_load", i_val=ival)
374 CALL dbcsr_set_config(comm_thread_load=ival)
375
376 CALL section_vals_val_get(dbcsr_section, "multrec_limit", i_val=ival)
377 CALL dbcsr_set_config(multrec_limit=ival)
378
379 CALL section_vals_val_get(dbcsr_section, "use_mempools_cpu", l_val=lval)
380 CALL dbcsr_set_config(use_mempools_cpu=lval)
381
382 CALL section_vals_val_get(dbcsr_section, "TENSOR%tas_split_factor", r_val=rval)
383 CALL dbcsr_set_config(tas_split_factor=rval)
384
385 CALL section_vals_val_get(dbcsr_section, "ACC%thread_buffers", i_val=ival)
386 CALL dbcsr_set_config(accdrv_thread_buffers=ival)
387
388 CALL section_vals_val_get(dbcsr_section, "ACC%min_flop_process", i_val=ival)
389 CALL dbcsr_set_config(accdrv_min_flop_process=ival)
390
391 CALL section_vals_val_get(dbcsr_section, "ACC%stack_sort", l_val=lval)
392 CALL dbcsr_set_config(accdrv_stack_sort=lval)
393
394 CALL section_vals_val_get(dbcsr_section, "ACC%min_flop_sort", i_val=ival)
395 CALL dbcsr_set_config(accdrv_min_flop_sort=ival)
396
397 CALL section_vals_val_get(dbcsr_section, "ACC%process_inhomogenous", l_val=lval)
398 CALL dbcsr_set_config(accdrv_do_inhomogenous=lval)
399
400 CALL section_vals_val_get(dbcsr_section, "ACC%avoid_after_busy", l_val=lval)
401 CALL dbcsr_set_config(accdrv_avoid_after_busy=lval)
402
403 CALL section_vals_val_get(dbcsr_section, "ACC%binning_nbins", i_val=ival)
404 CALL dbcsr_set_config(accdrv_binning_nbins=ival)
405
406 CALL section_vals_val_get(dbcsr_section, "ACC%binning_binsize", i_val=ival)
407 CALL dbcsr_set_config(accdrv_binning_binsize=ival)
408
409 CALL section_vals_val_get(dbcsr_section, "mm_driver", i_val=ival)
410 SELECT CASE (ival)
411 CASE (mm_driver_auto)
412 CALL dbcsr_set_config(mm_driver="AUTO")
413#if defined(__LIBXSMM)
414 CALL cite_reference(heinecke2016)
415#endif
416 CASE (mm_driver_blas)
417 CALL dbcsr_set_config(mm_driver="BLAS")
418 CASE (mm_driver_matmul)
419 CALL dbcsr_set_config(mm_driver="MATMUL")
420 CASE (mm_driver_smm)
421 CALL dbcsr_set_config(mm_driver="SMM")
422 CASE (mm_driver_xsmm)
423 CALL dbcsr_set_config(mm_driver="XSMM")
424 CALL cite_reference(heinecke2016)
425 CASE DEFAULT
426 cpabort("Unknown mm_driver")
427 END SELECT
428
429 CALL timestop(handle)
430 END SUBROUTINE cp_dbcsr_config
431
432! **************************************************************************************************
433!> \brief allocate the blocks of a dbcsr based on the neighbor list
434!> \param matrix the matrix
435!> \param sab_orb the corresponding neighbor list
436!> \param desymmetrize Allocate all block of a non-symmetric matrix from a symmetric list
437!> \par History
438!> 11.2009 created vw
439!> 01.2014 moved here from cp_dbcsr_operations (Ole Schuett)
440!> \author vw
441!> \note
442! **************************************************************************************************
443
444 SUBROUTINE cp_dbcsr_alloc_block_from_nbl(matrix, sab_orb, desymmetrize)
445
446 TYPE(dbcsr_type) :: matrix
447 TYPE(neighbor_list_set_p_type), DIMENSION(:), &
448 POINTER :: sab_orb
449 LOGICAL, INTENT(IN), OPTIONAL :: desymmetrize
450
451 CHARACTER(LEN=*), PARAMETER :: routinen = 'cp_dbcsr_alloc_block_from_nbl'
452
453 CHARACTER(LEN=1) :: symmetry
454 INTEGER :: blk_cnt, handle, iatom, icol, inode, &
455 irow, jatom, last_jatom, nadd
456 INTEGER, ALLOCATABLE, DIMENSION(:) :: cols, rows, tmp
457 LOGICAL :: alloc_full, is_symmetric, new_atom_b
459 DIMENSION(:), POINTER :: nl_iterator
460
461 CALL timeset(routinen, handle)
462
463 symmetry = dbcsr_get_matrix_type(matrix)
464
465 cpassert(ASSOCIATED(sab_orb))
466
467 CALL get_neighbor_list_set_p(neighbor_list_sets=sab_orb, symmetric=is_symmetric)
468 alloc_full = .false.
469 IF (PRESENT(desymmetrize)) THEN
470 IF (desymmetrize .AND. (symmetry == dbcsr_type_no_symmetry)) THEN
471 IF (is_symmetric) alloc_full = .true.
472 END IF
473 END IF
474
475 CALL dbcsr_finalize(matrix)
476 ALLOCATE (rows(3), cols(3))
477 blk_cnt = 0
478 nadd = 1
479 IF (alloc_full) nadd = 2
480
481 CALL neighbor_list_iterator_create(nl_iterator, sab_orb)
482 DO WHILE (neighbor_list_iterate(nl_iterator) == 0)
483 CALL get_iterator_info(nl_iterator, iatom=iatom, jatom=jatom, inode=inode)
484 IF (inode == 1) last_jatom = 0
485 IF (jatom /= last_jatom) THEN
486 new_atom_b = .true.
487 last_jatom = jatom
488 ELSE
489 new_atom_b = .false.
490 cycle
491 END IF
492 IF (blk_cnt + nadd .GT. SIZE(rows)) THEN
493 ALLOCATE (tmp(blk_cnt + nadd))
494 tmp(1:blk_cnt) = rows(1:blk_cnt)
495 DEALLOCATE (rows)
496 ALLOCATE (rows((blk_cnt + nadd)*2))
497 rows(1:blk_cnt) = tmp(1:blk_cnt)
498 tmp(1:blk_cnt) = cols(1:blk_cnt)
499 DEALLOCATE (cols)
500 ALLOCATE (cols((blk_cnt + nadd)*2))
501 cols(1:blk_cnt) = tmp(1:blk_cnt)
502 DEALLOCATE (tmp)
503 END IF
504 IF (alloc_full) THEN
505 blk_cnt = blk_cnt + 1
506 rows(blk_cnt) = iatom
507 cols(blk_cnt) = jatom
508 IF (iatom /= jatom) THEN
509 blk_cnt = blk_cnt + 1
510 rows(blk_cnt) = jatom
511 cols(blk_cnt) = iatom
512 END IF
513 ELSE
514 blk_cnt = blk_cnt + 1
515 IF (symmetry == dbcsr_type_no_symmetry) THEN
516 rows(blk_cnt) = iatom
517 cols(blk_cnt) = jatom
518 ELSE
519 IF (iatom <= jatom) THEN
520 irow = iatom
521 icol = jatom
522 ELSE
523 irow = jatom
524 icol = iatom
525 END IF
526 rows(blk_cnt) = irow
527 cols(blk_cnt) = icol
528 END IF
529 END IF
530
531 END DO
532 CALL neighbor_list_iterator_release(nl_iterator)
533
534 !
535 CALL dbcsr_reserve_blocks(matrix, rows(1:blk_cnt), cols(1:blk_cnt))
536 DEALLOCATE (rows)
537 DEALLOCATE (cols)
538 CALL dbcsr_finalize(matrix)
539
540 CALL timestop(handle)
541
542 END SUBROUTINE cp_dbcsr_alloc_block_from_nbl
543
544! **************************************************************************************************
545!> \brief Apply distance screening to refine sparsity pattern of matrices in CSR
546!> format (using eps_pgf_orb). Currently this is used for the external
547!> library PEXSI.
548!> \param ks_env ...
549!> \param[in, out] csr_sparsity DBCSR matrix defining CSR sparsity pattern.
550!> This matrix must be initialized and allocated
551!> with exactly the same DBCSR sparsity pattern as
552!> the DBCSR matrix that is used to create the CSR
553!> matrix. It must have symmetric DBCSR format and
554!> must not be filtered.
555!> \par History
556!> 02.2015 created [Patrick Seewald]
557!> \author Patrick Seewald
558! **************************************************************************************************
559 SUBROUTINE cp_dbcsr_to_csr_screening(ks_env, csr_sparsity)
560 TYPE(qs_ks_env_type), POINTER :: ks_env
561 TYPE(dbcsr_type), INTENT(INOUT) :: csr_sparsity
562
563 CHARACTER(len=*), PARAMETER :: routinen = 'cp_dbcsr_to_csr_screening'
564
565 INTEGER :: atom_a, atom_b, handle, iatom, icol, ikind, ipgf, irow, iset, isgf, ishell, &
566 jatom, jkind, jpgf, jset, jsgf, jshell, nkind, nset_a, nset_b
567 INTEGER, ALLOCATABLE, DIMENSION(:) :: atom_of_kind
568 INTEGER, DIMENSION(:), POINTER :: npgf_a, npgf_b, nshell_a, nshell_b
569 INTEGER, DIMENSION(:, :), POINTER :: l_a, l_b
570 LOGICAL :: do_symmetric, found
571 REAL(kind=dp) :: dab, eps_pgf_orb, r_a, r_b
572 REAL(kind=dp), DIMENSION(3) :: rab
573 REAL(kind=dp), DIMENSION(:, :), POINTER :: rpgfa, rpgfb, zet_a, zet_b
574 REAL(kind=dp), DIMENSION(:, :, :), POINTER :: gcc_a, gcc_b
575 REAL(kind=real_8), DIMENSION(:, :), POINTER :: screen_blk
576 TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
577 TYPE(dft_control_type), POINTER :: dft_control
578 TYPE(gto_basis_set_p_type), DIMENSION(:), POINTER :: basis_set_list_a, basis_set_list_b
579 TYPE(gto_basis_set_type), POINTER :: basis_set_a, basis_set_b
581 DIMENSION(:), POINTER :: nl_iterator
582 TYPE(neighbor_list_set_p_type), DIMENSION(:), &
583 POINTER :: neighbour_list
584 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
585
586 NULLIFY (screen_blk, atomic_kind_set, basis_set_list_a, &
587 basis_set_list_b, basis_set_a, basis_set_b, nl_iterator, &
588 qs_kind_set, dft_control)
589
590 CALL timeset(routinen, handle)
591
592 cpassert(dbcsr_has_symmetry(csr_sparsity))
593
594 CALL get_ks_env(ks_env, &
595 sab_orb=neighbour_list, &
596 atomic_kind_set=atomic_kind_set, &
597 qs_kind_set=qs_kind_set, &
598 dft_control=dft_control)
599
600 eps_pgf_orb = dft_control%qs_control%eps_pgf_orb
601
602 nkind = SIZE(qs_kind_set)
603 CALL get_atomic_kind_set(atomic_kind_set, atom_of_kind=atom_of_kind)
604 cpassert(SIZE(neighbour_list) > 0)
605 CALL get_neighbor_list_set_p(neighbor_list_sets=neighbour_list, symmetric=do_symmetric)
606 cpassert(do_symmetric)
607 ALLOCATE (basis_set_list_a(nkind), basis_set_list_b(nkind))
608 CALL basis_set_list_setup(basis_set_list_a, "ORB", qs_kind_set)
609 CALL basis_set_list_setup(basis_set_list_b, "ORB", qs_kind_set)
610
611 ! csr_sparsity can obtain values 0 (if zero element) or 1 (if non-zero element)
612 CALL dbcsr_set(csr_sparsity, 0.0)
613
614 CALL neighbor_list_iterator_create(nl_iterator, neighbour_list)
615
616 ! Iterate over interacting pairs of atoms corresponding to non-zero
617 ! DBCSR blocks
618 DO WHILE (neighbor_list_iterate(nl_iterator) == 0)
619 CALL get_iterator_info(nl_iterator, &
620 ikind=ikind, jkind=jkind, &
621 iatom=iatom, jatom=jatom, &
622 r=rab)
623
624 basis_set_a => basis_set_list_a(ikind)%gto_basis_set
625 IF (.NOT. ASSOCIATED(basis_set_a)) cycle
626 basis_set_b => basis_set_list_b(jkind)%gto_basis_set
627 IF (.NOT. ASSOCIATED(basis_set_b)) cycle
628
629 atom_a = atom_of_kind(iatom)
630 atom_b = atom_of_kind(jatom)
631
632 nset_a = basis_set_a%nset
633 nset_b = basis_set_b%nset
634 npgf_a => basis_set_a%npgf
635 npgf_b => basis_set_b%npgf
636 nshell_a => basis_set_a%nshell
637 nshell_b => basis_set_b%nshell
638
639 l_a => basis_set_a%l
640 l_b => basis_set_b%l
641 gcc_a => basis_set_a%gcc
642 gcc_b => basis_set_b%gcc
643 zet_a => basis_set_a%zet
644 zet_b => basis_set_b%zet
645
646 rpgfa => basis_set_a%pgf_radius
647 rpgfb => basis_set_b%pgf_radius
648
649 IF (iatom <= jatom) THEN
650 irow = iatom
651 icol = jatom
652 ELSE
653 irow = jatom
654 icol = iatom
655 END IF
656
657 CALL dbcsr_get_block_p(matrix=csr_sparsity, row=irow, col=icol, &
658 block=screen_blk, found=found)
659
660 cpassert(found)
661
662 ! Distance between atoms a and b
663 dab = sqrt(rab(1)**2 + rab(2)**2 + rab(3)**2)
664
665 ! iterate over pairs of primitive GTOs i,j, get their radii r_i, r_j according
666 ! to eps_pgf_orb. Define all matrix elements as non-zero to which a
667 ! contribution from two Gaussians i,j exists with r_i + r_j >= dab.
668
669 isgf = 0
670 DO iset = 1, nset_a
671 DO ishell = 1, nshell_a(iset)
672 jsgf = 0
673 DO jset = 1, nset_b
674 DO jshell = 1, nshell_b(jset)
675 gto_loop: DO ipgf = 1, npgf_a(iset)
676 DO jpgf = 1, npgf_b(jset)
677 IF (rpgfa(ipgf, iset) + rpgfb(jpgf, jset) .GE. dab) THEN
678 ! more selective screening with radius calculated for each primitive GTO
679 r_a = exp_radius(l_a(ishell, iset), &
680 zet_a(ipgf, iset), &
681 eps_pgf_orb, &
682 gcc_a(ipgf, ishell, iset))
683 r_b = exp_radius(l_b(jshell, jset), &
684 zet_b(jpgf, jset), &
685 eps_pgf_orb, &
686 gcc_b(jpgf, jshell, jset))
687 IF (r_a + r_b .GE. dab) THEN
688 IF (irow .EQ. iatom) THEN
689 screen_blk(isgf + 1:isgf + nso(l_a(ishell, iset)), &
690 jsgf + 1:jsgf + nso(l_b(jshell, jset))) = 1.0_dp
691 ELSE
692 screen_blk(jsgf + 1:jsgf + nso(l_b(jshell, jset)), &
693 isgf + 1:isgf + nso(l_a(ishell, iset))) = 1.0_dp
694 END IF
695 EXIT gto_loop
696 END IF
697 END IF
698 END DO
699 END DO gto_loop
700 jsgf = jsgf + nso(l_b(jshell, jset))
701 END DO
702 END DO
703 isgf = isgf + nso(l_a(ishell, iset))
704 END DO
705 END DO
706 END DO
707
708 CALL neighbor_list_iterator_release(nl_iterator)
709 DEALLOCATE (basis_set_list_a, basis_set_list_b)
710
711 CALL timestop(handle)
712 END SUBROUTINE cp_dbcsr_to_csr_screening
713
714END MODULE cp_dbcsr_cp2k_link
All kind of helpful little routines.
Definition ao_util.F:14
real(kind=dp) function, public exp_radius(l, alpha, threshold, prefactor, epsabs, epsrel, rlow)
The radius of a primitive Gaussian function for a given threshold is calculated. g(r) = prefactor*r**...
Definition ao_util.F:96
Define the atomic kind types and their sub types.
subroutine, public get_atomic_kind_set(atomic_kind_set, atom_of_kind, kind_of, natom_of_kind, maxatom, natom, nshell, fist_potential_present, shell_present, shell_adiabatic, shell_check_distance, damping_present)
Get attributes of an atomic kind set.
collects all references to literature in CP2K as new algorithms / method are included from literature...
integer, save, public heinecke2016
integer, save, public borstnik2014
integer, save, public schuett2016
Defines control structures, which contain the parameters and the settings for the DFT-based calculati...
DBCSR operations in CP2K.
integer, save, public max_elements_per_block
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 type(section_vals_type) function, pointer, public section_vals_get_subs_vals(section_vals, subsection_name, i_rep_section, can_return_null)
returns the values of the requested subsection
recursive subroutine, public section_release(section)
releases the given keyword list (see doc/ReferenceCounting.html)
subroutine, public section_vals_val_get(section_vals, keyword_name, i_rep_section, i_rep_val, n_rep_val, val, l_val, i_val, r_val, c_val, l_vals, i_vals, r_vals, c_vals, explicit)
returns the requested value
Defines the basic variable types.
Definition kinds.F:23
integer, parameter, public dp
Definition kinds.F:34
integer, parameter, public real_8
Definition kinds.F:41
Provides Cartesian and spherical orbital pointers and indices.
integer, dimension(:), allocatable, public nso
Some utility functions for the calculation of integrals.
subroutine, public basis_set_list_setup(basis_set_list, basis_type, qs_kind_set)
Set up an easy accessible list of the basis sets for all kinds.
Define the quickstep kind type and their sub types.
subroutine, public get_ks_env(ks_env, v_hartree_rspace, s_mstruct_changed, rho_changed, potential_changed, forces_up_to_date, complex_ks, matrix_h, matrix_h_im, matrix_ks, matrix_ks_im, matrix_vxc, kinetic, matrix_s, matrix_s_ri_aux, matrix_w, matrix_p_mp2, matrix_p_mp2_admm, matrix_h_kp, matrix_h_im_kp, matrix_ks_kp, matrix_vxc_kp, kinetic_kp, matrix_s_kp, matrix_w_kp, matrix_s_ri_aux_kp, matrix_ks_im_kp, rho, rho_xc, vppl, rho_core, rho_nlcc, rho_nlcc_g, vee, neighbor_list_id, sab_orb, sab_all, sac_ae, sac_ppl, sac_lri, sap_ppnl, sap_oce, sab_lrc, sab_se, sab_xtbe, sab_tbe, sab_core, sab_xb, sab_xtb_nonbond, sab_vdw, sab_scp, sab_almo, sab_kp, sab_kp_nosym, task_list, task_list_soft, kpoints, do_kpoints, atomic_kind_set, qs_kind_set, cell, cell_ref, use_ref_cell, particle_set, energy, force, local_particles, local_molecules, molecule_kind_set, molecule_set, subsys, cp_subsys, virial, results, atprop, nkind, natom, dft_control, dbcsr_dist, distribution_2d, pw_env, para_env, blacs_env, nelectron_total, nelectron_spin)
...
Define the neighbor list data types and the corresponding functionality.
subroutine, public neighbor_list_iterator_create(iterator_set, nl, search, nthread)
Neighbor list iterator functions.
subroutine, public neighbor_list_iterator_release(iterator_set)
...
subroutine, public get_neighbor_list_set_p(neighbor_list_sets, nlist, symmetric)
Return the components of the first neighbor list set.
integer function, public neighbor_list_iterate(iterator_set, mepos)
...
subroutine, public get_iterator_info(iterator_set, mepos, ikind, jkind, nkind, ilist, nlist, inode, nnode, iatom, jatom, r, cell)
...
Utilities for string manipulations.
Provides all information about an atomic kind.
represent a keyword in the input
represent a section of the input file
Provides all information about a quickstep kind.
calculation environment to calculate the ks matrix, holds all the needed vars. assumes that the core ...