(git:d18deda)
Loading...
Searching...
No Matches
cp_dbcsr_operations.F
Go to the documentation of this file.
1!--------------------------------------------------------------------------------------------------!
2! CP2K: A general program to perform molecular dynamics simulations !
3! Copyright 2000-2025 CP2K developers group <https://cp2k.org> !
4! !
5! SPDX-License-Identifier: GPL-2.0-or-later !
6!--------------------------------------------------------------------------------------------------!
7
8! **************************************************************************************************
9!> \brief DBCSR operations in CP2K
10!> \author Urban Borstnik
11!> \date 2009-05-12
12!> \version 0.8
13!>
14!> <b>Modification history:</b>
15!> - Created 2009-05-12
16!> - Generalized sm_fm_mulitply for matrices w/ different row/col block size (A. Bussy, 11.2018)
17! **************************************************************************************************
18MODULE cp_dbcsr_operations
19 USE cp_blacs_env, ONLY: cp_blacs_env_type
20 USE cp_dbcsr_api, ONLY: &
21 dbcsr_add, dbcsr_complete_redistribute, dbcsr_convert_sizes_to_offsets, dbcsr_copy, &
22 dbcsr_create, dbcsr_deallocate_matrix, dbcsr_desymmetrize, dbcsr_distribution_get, &
23 dbcsr_distribution_new, dbcsr_distribution_release, dbcsr_distribution_type, &
24 dbcsr_get_info, dbcsr_get_matrix_type, dbcsr_iterator_blocks_left, &
25 dbcsr_iterator_next_block, dbcsr_iterator_readonly_start, dbcsr_iterator_start, &
26 dbcsr_iterator_stop, dbcsr_iterator_type, dbcsr_multiply, dbcsr_p_type, dbcsr_release, &
27 dbcsr_scale, dbcsr_type, dbcsr_type_antisymmetric, dbcsr_type_no_symmetry, &
28 dbcsr_type_symmetric, dbcsr_valid_index, dbcsr_verify_matrix
29 USE cp_dbcsr_contrib, ONLY: dbcsr_frobenius_norm,&
30 dbcsr_reserve_all_blocks
31 USE cp_fm_basic_linalg, ONLY: cp_fm_gemm
32 USE cp_fm_struct, ONLY: cp_fm_struct_create,&
33 cp_fm_struct_release,&
34 cp_fm_struct_type
35 USE cp_fm_types, ONLY: cp_fm_create,&
36 cp_fm_get_info,&
37 cp_fm_release,&
38 cp_fm_to_fm,&
39 cp_fm_type
40 USE distribution_2d_types, ONLY: distribution_2d_get,&
41 distribution_2d_type
42 USE kinds, ONLY: default_string_length,&
43 dp
44 USE mathlib, ONLY: gcd,&
45 lcm
46 USE message_passing, ONLY: mp_para_env_type
47
48!$ USE OMP_LIB, ONLY: omp_get_max_threads, omp_get_thread_num, omp_get_num_threads
49#include "base/base_uses.f90"
50
51 IMPLICIT NONE
52 PRIVATE
53
54 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'cp_dbcsr_operations'
55 LOGICAL, PARAMETER :: debug_mod = .false.
56
57 INTEGER, SAVE, PUBLIC :: max_elements_per_block = 32
58
59 PUBLIC :: dbcsr_multiply_local
60
61 ! CP2K API emulation
62 PUBLIC :: copy_fm_to_dbcsr, copy_dbcsr_to_fm, &
63 cp_dbcsr_sm_fm_multiply, cp_dbcsr_plus_fm_fm_t, &
64 copy_dbcsr_to_fm_bc, copy_fm_to_dbcsr_bc, cp_fm_to_dbcsr_row_template, &
65 cp_dbcsr_m_by_n_from_template, cp_dbcsr_m_by_n_from_row_template, &
66 dbcsr_create_dist_r_unrot
67
68 ! distribution_2d_type compatibility
69 PUBLIC :: cp_dbcsr_dist2d_to_dist
70
71 PUBLIC :: dbcsr_copy_columns_hack
72
73 ! matrix set
74 PUBLIC :: dbcsr_allocate_matrix_set
75 PUBLIC :: dbcsr_deallocate_matrix_set
76
77 INTERFACE dbcsr_allocate_matrix_set
78 MODULE PROCEDURE allocate_dbcsr_matrix_set_1d
79 MODULE PROCEDURE allocate_dbcsr_matrix_set_2d
80 MODULE PROCEDURE allocate_dbcsr_matrix_set_3d
81 MODULE PROCEDURE allocate_dbcsr_matrix_set_4d
82 MODULE PROCEDURE allocate_dbcsr_matrix_set_5d
83 END INTERFACE
84
85 INTERFACE dbcsr_deallocate_matrix_set
86 MODULE PROCEDURE deallocate_dbcsr_matrix_set_1d
87 MODULE PROCEDURE deallocate_dbcsr_matrix_set_2d
88 MODULE PROCEDURE deallocate_dbcsr_matrix_set_3d
89 MODULE PROCEDURE deallocate_dbcsr_matrix_set_4d
90 MODULE PROCEDURE deallocate_dbcsr_matrix_set_5d
91 END INTERFACE
92
93CONTAINS
94
95! **************************************************************************************************
96!> \brief Copy a BLACS matrix to a dbcsr matrix.
97!>
98!> real_matrix=beta*real_matrix+alpha*fm
99!> beta defaults to 0, alpha to 1
100!> \param[in] fm full matrix
101!> \param[out] matrix DBCSR matrix
102!> \param[in] keep_sparsity (optional) retains the sparsity of the input
103!> matrix
104!> \date 2009-10-13
105!> \par History
106!> 2009-10-13 rewritten based on copy_dbcsr_to_fm
107!> \author Urban Borstnik
108!> \version 2.0
109! **************************************************************************************************
110 SUBROUTINE copy_fm_to_dbcsr(fm, matrix, keep_sparsity)
111 TYPE(cp_fm_type), INTENT(IN) :: fm
112 TYPE(dbcsr_type), INTENT(INOUT) :: matrix
113 LOGICAL, INTENT(IN), OPTIONAL :: keep_sparsity
114
115 CHARACTER(LEN=*), PARAMETER :: routinen = 'copy_fm_to_dbcsr'
116
117 INTEGER :: handle
118 LOGICAL :: my_keep_sparsity
119 TYPE(dbcsr_type) :: bc_mat, redist_mat
120
121 CALL timeset(routinen, handle)
122
123 my_keep_sparsity = .false.
124 IF (PRESENT(keep_sparsity)) my_keep_sparsity = keep_sparsity
125
126 CALL copy_fm_to_dbcsr_bc(fm, bc_mat)
127
128 IF (my_keep_sparsity) THEN
129 CALL dbcsr_create(redist_mat, template=matrix)
130 CALL dbcsr_complete_redistribute(bc_mat, redist_mat)
131 CALL dbcsr_copy(matrix, redist_mat, keep_sparsity=.true.)
132 CALL dbcsr_release(redist_mat)
133 ELSE
134 CALL dbcsr_complete_redistribute(bc_mat, matrix)
135 END IF
136
137 CALL dbcsr_release(bc_mat)
138
139 CALL timestop(handle)
140 END SUBROUTINE copy_fm_to_dbcsr
141
142! **************************************************************************************************
143!> \brief Copy a BLACS matrix to a dbcsr matrix with a special block-cyclic distribution,
144!> which requires no complete redistribution.
145!> \param fm ...
146!> \param bc_mat ...
147! **************************************************************************************************
148 SUBROUTINE copy_fm_to_dbcsr_bc(fm, bc_mat)
149 TYPE(cp_fm_type), INTENT(IN) :: fm
150 TYPE(dbcsr_type), INTENT(INOUT) :: bc_mat
151
152 CHARACTER(LEN=*), PARAMETER :: routinen = 'copy_fm_to_dbcsr_bc'
153
154 INTEGER :: col, handle, ncol_block, ncol_global, &
155 nrow_block, nrow_global, row
156 INTEGER, ALLOCATABLE, DIMENSION(:) :: first_col, first_row, last_col, last_row
157 INTEGER, DIMENSION(:), POINTER :: col_blk_size, row_blk_size
158 INTEGER, DIMENSION(:, :), POINTER :: pgrid
159 REAL(kind=dp), DIMENSION(:, :), POINTER :: dbcsr_block, fm_block
160 TYPE(dbcsr_distribution_type) :: bc_dist
161 TYPE(dbcsr_iterator_type) :: iter
162
163 CALL timeset(routinen, handle)
164
165 IF (fm%use_sp) cpabort("copy_fm_to_dbcsr_bc: single precision not supported")
166
167 ! Create processor grid
168 pgrid => fm%matrix_struct%context%blacs2mpi
169
170 ! Create a block-cyclic distribution compatible with the FM matrix.
171 nrow_block = fm%matrix_struct%nrow_block
172 ncol_block = fm%matrix_struct%ncol_block
173 nrow_global = fm%matrix_struct%nrow_global
174 ncol_global = fm%matrix_struct%ncol_global
175 NULLIFY (col_blk_size, row_blk_size)
176 CALL dbcsr_create_dist_block_cyclic(bc_dist, &
177 nrows=nrow_global, ncolumns=ncol_global, & ! Actual full matrix size
178 nrow_block=nrow_block, ncol_block=ncol_block, & ! BLACS parameters
179 group_handle=fm%matrix_struct%para_env%get_handle(), pgrid=pgrid, &
180 row_blk_sizes=row_blk_size, col_blk_sizes=col_blk_size) ! block-cyclic row/col sizes
181
182 ! Create the block-cyclic DBCSR matrix
183 CALL dbcsr_create(bc_mat, "Block-cyclic ", bc_dist, &
184 dbcsr_type_no_symmetry, row_blk_size, col_blk_size, reuse_arrays=.true.)
185 CALL dbcsr_distribution_release(bc_dist)
186
187 ! allocate all blocks
188 CALL dbcsr_reserve_all_blocks(bc_mat)
189
190 CALL calculate_fm_block_ranges(bc_mat, first_row, last_row, first_col, last_col)
191
192 ! Copy the FM data to the block-cyclic DBCSR matrix. This step
193 ! could be skipped with appropriate DBCSR index manipulation.
194 fm_block => fm%local_data
195!$OMP PARALLEL DEFAULT(NONE) PRIVATE(iter, row, col, dbcsr_block) &
196!$OMP SHARED(bc_mat, last_row, first_row, last_col, first_col, fm_block)
197 CALL dbcsr_iterator_start(iter, bc_mat)
198 DO WHILE (dbcsr_iterator_blocks_left(iter))
199 CALL dbcsr_iterator_next_block(iter, row, col, dbcsr_block)
200 dbcsr_block(:, :) = fm_block(first_row(row):last_row(row), first_col(col):last_col(col))
201 END DO
202 CALL dbcsr_iterator_stop(iter)
203!$OMP END PARALLEL
204
205 CALL timestop(handle)
206 END SUBROUTINE copy_fm_to_dbcsr_bc
207
208! **************************************************************************************************
209!> \brief Copy a DBCSR matrix to a BLACS matrix
210!> \param[in] matrix DBCSR matrix
211!> \param[out] fm full matrix
212! **************************************************************************************************
213 SUBROUTINE copy_dbcsr_to_fm(matrix, fm)
214 TYPE(dbcsr_type), INTENT(IN) :: matrix
215 TYPE(cp_fm_type), INTENT(INOUT) :: fm
216
217 CHARACTER(LEN=*), PARAMETER :: routinen = 'copy_dbcsr_to_fm'
218
219 CHARACTER(len=default_string_length) :: name
220 INTEGER :: group_handle, handle, ncol_block, &
221 nfullcols_total, nfullrows_total, &
222 nrow_block
223 INTEGER, DIMENSION(:), POINTER :: col_blk_size, row_blk_size
224 INTEGER, DIMENSION(:, :), POINTER :: pgrid
225 TYPE(dbcsr_distribution_type) :: bc_dist, dist
226 TYPE(dbcsr_type) :: bc_mat, matrix_nosym
227
228 CALL timeset(routinen, handle)
229
230 ! check compatibility
231 CALL dbcsr_get_info(matrix, &
232 name=name, &
233 distribution=dist, &
234 nfullrows_total=nfullrows_total, &
235 nfullcols_total=nfullcols_total)
236
237 cpassert(fm%matrix_struct%nrow_global == nfullrows_total)
238 cpassert(fm%matrix_struct%ncol_global == nfullcols_total)
239
240 ! info about the full matrix
241 nrow_block = fm%matrix_struct%nrow_block
242 ncol_block = fm%matrix_struct%ncol_block
243
244 ! Convert DBCSR to a block-cyclic
245 NULLIFY (col_blk_size, row_blk_size)
246 CALL dbcsr_distribution_get(dist, group=group_handle, pgrid=pgrid)
247 CALL dbcsr_create_dist_block_cyclic(bc_dist, &
248 nrows=nfullrows_total, ncolumns=nfullcols_total, &
249 nrow_block=nrow_block, ncol_block=ncol_block, &
250 group_handle=group_handle, pgrid=pgrid, &
251 row_blk_sizes=row_blk_size, col_blk_sizes=col_blk_size)
252
253 CALL dbcsr_create(bc_mat, "Block-cyclic"//name, bc_dist, &
254 dbcsr_type_no_symmetry, row_blk_size, col_blk_size, reuse_arrays=.true.)
255 CALL dbcsr_distribution_release(bc_dist)
256
257 CALL dbcsr_create(matrix_nosym, template=matrix, matrix_type="N")
258 CALL dbcsr_desymmetrize(matrix, matrix_nosym)
259 CALL dbcsr_complete_redistribute(matrix_nosym, bc_mat)
260 CALL dbcsr_release(matrix_nosym)
261
262 CALL copy_dbcsr_to_fm_bc(bc_mat, fm)
263
264 CALL dbcsr_release(bc_mat)
265
266 CALL timestop(handle)
267 END SUBROUTINE copy_dbcsr_to_fm
268
269! **************************************************************************************************
270!> \brief Copy a DBCSR_BLACS matrix to a BLACS matrix
271!> \param bc_mat DBCSR matrix
272!> \param[out] fm full matrix
273! **************************************************************************************************
274 SUBROUTINE copy_dbcsr_to_fm_bc(bc_mat, fm)
275 TYPE(dbcsr_type), INTENT(IN) :: bc_mat
276 TYPE(cp_fm_type), INTENT(INOUT) :: fm
277
278 CHARACTER(LEN=*), PARAMETER :: routinen = 'copy_dbcsr_to_fm_bc'
279
280 INTEGER :: col, handle, row
281 INTEGER, ALLOCATABLE, DIMENSION(:) :: first_col, first_row, last_col, last_row
282 REAL(kind=dp), DIMENSION(:, :), POINTER :: dbcsr_block, fm_block
283 TYPE(dbcsr_iterator_type) :: iter
284
285 CALL timeset(routinen, handle)
286
287 IF (fm%use_sp) cpabort("copy_dbcsr_to_fm_bc: single precision not supported")
288
289 CALL calculate_fm_block_ranges(bc_mat, first_row, last_row, first_col, last_col)
290
291 ! Now copy data to the FM matrix
292 fm_block => fm%local_data
293 fm_block = real(0.0, kind=dp)
294!$OMP PARALLEL DEFAULT(NONE) PRIVATE(iter, row, col, dbcsr_block) &
295!$OMP SHARED(bc_mat, last_row, first_row, last_col, first_col, fm_block)
296 CALL dbcsr_iterator_readonly_start(iter, bc_mat)
297 DO WHILE (dbcsr_iterator_blocks_left(iter))
298 CALL dbcsr_iterator_next_block(iter, row, col, dbcsr_block)
299 fm_block(first_row(row):last_row(row), first_col(col):last_col(col)) = dbcsr_block(:, :)
300 END DO
301 CALL dbcsr_iterator_stop(iter)
302!$OMP END PARALLEL
303
304 CALL timestop(handle)
305 END SUBROUTINE copy_dbcsr_to_fm_bc
306
307! **************************************************************************************************
308!> \brief Helper routine used to copy blocks from DBCSR into FM matrices and vice versa
309!> \param bc_mat ...
310!> \param first_row ...
311!> \param last_row ...
312!> \param first_col ...
313!> \param last_col ...
314!> \author Ole Schuett
315! **************************************************************************************************
316 SUBROUTINE calculate_fm_block_ranges(bc_mat, first_row, last_row, first_col, last_col)
317 TYPE(dbcsr_type), INTENT(IN) :: bc_mat
318 INTEGER, ALLOCATABLE, DIMENSION(:), INTENT(OUT) :: first_row, last_row, first_col, last_col
319
320 INTEGER :: col, nblkcols_local, nblkcols_total, &
321 nblkrows_local, nblkrows_total, row
322 INTEGER, ALLOCATABLE, DIMENSION(:) :: local_col_sizes, local_row_sizes
323 INTEGER, DIMENSION(:), POINTER :: col_blk_size, local_cols, local_rows, &
324 row_blk_size
325
326 CALL dbcsr_get_info(bc_mat, &
327 nblkrows_total=nblkrows_total, &
328 nblkcols_total=nblkcols_total, &
329 nblkrows_local=nblkrows_local, &
330 nblkcols_local=nblkcols_local, &
331 local_rows=local_rows, &
332 local_cols=local_cols, &
333 row_blk_size=row_blk_size, &
334 col_blk_size=col_blk_size)
335
336 ! calculate first_row and last_row
337 ALLOCATE (local_row_sizes(nblkrows_total))
338 local_row_sizes(:) = 0
339 IF (nblkrows_local .GE. 1) THEN
340 DO row = 1, nblkrows_local
341 local_row_sizes(local_rows(row)) = row_blk_size(local_rows(row))
342 END DO
343 END IF
344 ALLOCATE (first_row(nblkrows_total), last_row(nblkrows_total))
345 CALL dbcsr_convert_sizes_to_offsets(local_row_sizes, first_row, last_row)
346 DEALLOCATE (local_row_sizes)
347
348 ! calculate first_col and last_col
349 ALLOCATE (local_col_sizes(nblkcols_total))
350 local_col_sizes(:) = 0
351 IF (nblkcols_local .GE. 1) THEN
352 DO col = 1, nblkcols_local
353 local_col_sizes(local_cols(col)) = col_blk_size(local_cols(col))
354 END DO
355 END IF
356 ALLOCATE (first_col(nblkcols_total), last_col(nblkcols_total))
357 CALL dbcsr_convert_sizes_to_offsets(local_col_sizes, first_col, last_col)
358 DEALLOCATE (local_col_sizes)
359
360 END SUBROUTINE calculate_fm_block_ranges
361
362! **************************************************************************************************
363!> \brief hack for dbcsr_copy_columns
364!> \param matrix_b ...
365!> \param matrix_a ...
366!> \param ncol ...
367!> \param source_start ...
368!> \param target_start ...
369!> \param para_env ...
370!> \param blacs_env ...
371!> \author vw
372! **************************************************************************************************
373 SUBROUTINE dbcsr_copy_columns_hack(matrix_b, matrix_a, &
374 ncol, source_start, target_start, para_env, blacs_env)
375
376 TYPE(dbcsr_type), INTENT(INOUT) :: matrix_b
377 TYPE(dbcsr_type), INTENT(IN) :: matrix_a
378 INTEGER, INTENT(IN) :: ncol, source_start, target_start
379 TYPE(mp_para_env_type), POINTER :: para_env
380 TYPE(cp_blacs_env_type), POINTER :: blacs_env
381
382 INTEGER :: nfullcols_total, nfullrows_total
383 TYPE(cp_fm_struct_type), POINTER :: fm_struct
384 TYPE(cp_fm_type) :: fm_matrix_a, fm_matrix_b
385
386 NULLIFY (fm_struct)
387 CALL dbcsr_get_info(matrix_a, nfullrows_total=nfullrows_total, nfullcols_total=nfullcols_total)
388 CALL cp_fm_struct_create(fm_struct, context=blacs_env, nrow_global=nfullrows_total, &
389 ncol_global=nfullcols_total, para_env=para_env)
390 CALL cp_fm_create(fm_matrix_a, fm_struct, name="fm_matrix_a")
391 CALL cp_fm_struct_release(fm_struct)
392
393 CALL dbcsr_get_info(matrix_b, nfullrows_total=nfullrows_total, nfullcols_total=nfullcols_total)
394 CALL cp_fm_struct_create(fm_struct, context=blacs_env, nrow_global=nfullrows_total, &
395 ncol_global=nfullcols_total, para_env=para_env)
396 CALL cp_fm_create(fm_matrix_b, fm_struct, name="fm_matrix_b")
397 CALL cp_fm_struct_release(fm_struct)
398
399 CALL copy_dbcsr_to_fm(matrix_a, fm_matrix_a)
400 CALL copy_dbcsr_to_fm(matrix_b, fm_matrix_b)
401
402 CALL cp_fm_to_fm(fm_matrix_a, fm_matrix_b, ncol, source_start, target_start)
403
404 CALL copy_fm_to_dbcsr(fm_matrix_b, matrix_b)
405
406 CALL cp_fm_release(fm_matrix_a)
407 CALL cp_fm_release(fm_matrix_b)
408
409 END SUBROUTINE dbcsr_copy_columns_hack
410
411! **************************************************************************************************
412!> \brief Creates a DBCSR distribution from a distribution_2d
413!> \param[in] dist2d distribution_2d
414!> \param[out] dist DBCSR distribution
415!> \par History
416!> move form dbcsr_operation 01.2010
417! **************************************************************************************************
418 SUBROUTINE cp_dbcsr_dist2d_to_dist(dist2d, dist)
419 TYPE(distribution_2d_type), INTENT(IN), TARGET :: dist2d
420 TYPE(dbcsr_distribution_type), INTENT(OUT) :: dist
421
422 INTEGER, DIMENSION(:), POINTER :: col_dist, row_dist
423 INTEGER, DIMENSION(:, :), POINTER :: col_dist_data, pgrid, row_dist_data
424 TYPE(cp_blacs_env_type), POINTER :: blacs_env
425 TYPE(distribution_2d_type), POINTER :: dist2d_p
426 TYPE(mp_para_env_type), POINTER :: para_env
427
428 dist2d_p => dist2d
429 CALL distribution_2d_get(dist2d_p, &
430 row_distribution=row_dist_data, &
431 col_distribution=col_dist_data, &
432 blacs_env=blacs_env)
433 CALL blacs_env%get(para_env=para_env, blacs2mpi=pgrid)
434
435 ! map to 1D arrays
436 row_dist => row_dist_data(:, 1)
437 col_dist => col_dist_data(:, 1)
438 !row_cluster => row_dist_data(:, 2)
439 !col_cluster => col_dist_data(:, 2)
440
441 CALL dbcsr_distribution_new(dist, &
442 group=para_env%get_handle(), pgrid=pgrid, &
443 row_dist=row_dist, &
444 col_dist=col_dist)
445
446 END SUBROUTINE cp_dbcsr_dist2d_to_dist
447
448! **************************************************************************************************
449!> \brief multiply a dbcsr with a replicated array
450!> c = alpha_scalar * A (dbscr) * b + c
451!> \param[in] matrix_a DBSCR matrxx
452!> \param[in] vec_b vectors b
453!> \param[inout] vec_c vectors c
454!> \param[in] ncol nbr of columns
455!> \param[in] alpha alpha
456!>
457! **************************************************************************************************
458 SUBROUTINE dbcsr_multiply_local(matrix_a, vec_b, vec_c, ncol, alpha)
459 TYPE(dbcsr_type), INTENT(IN) :: matrix_a
460 REAL(dp), DIMENSION(:, :), INTENT(IN) :: vec_b
461 REAL(dp), DIMENSION(:, :), INTENT(INOUT) :: vec_c
462 INTEGER, INTENT(in), OPTIONAL :: ncol
463 REAL(dp), INTENT(IN), OPTIONAL :: alpha
464
465 CHARACTER(LEN=*), PARAMETER :: routinen = 'dbcsr_multiply_local'
466
467 INTEGER :: col, coloff, my_ncol, row, rowoff, &
468 timing_handle
469 LOGICAL :: has_symm
470 REAL(dp) :: my_alpha, my_alpha2
471 REAL(dp), DIMENSION(:, :), POINTER :: data_d
472 TYPE(dbcsr_iterator_type) :: iter
473
474 CALL timeset(routinen, timing_handle)
475
476 my_alpha = 1.0_dp
477 IF (PRESENT(alpha)) my_alpha = alpha
478
479 my_ncol = SIZE(vec_b, 2)
480 IF (PRESENT(ncol)) my_ncol = ncol
481
482 my_alpha2 = 0.0_dp
483 IF (dbcsr_get_matrix_type(matrix_a) .EQ. dbcsr_type_symmetric) my_alpha2 = my_alpha
484 IF (dbcsr_get_matrix_type(matrix_a) .EQ. dbcsr_type_antisymmetric) my_alpha2 = -my_alpha
485
486 has_symm = (dbcsr_get_matrix_type(matrix_a) .EQ. dbcsr_type_symmetric .OR. &
487 dbcsr_get_matrix_type(matrix_a) .EQ. dbcsr_type_antisymmetric)
488
489!$OMP PARALLEL DEFAULT(NONE) SHARED(matrix_a,vec_b,vec_c,ncol,my_alpha2,my_alpha,my_ncol,has_symm) &
490!$OMP PRIVATE(iter,row,col,data_d,rowoff,coloff)
491 CALL dbcsr_iterator_readonly_start(iter, matrix_a, dynamic=.true., dynamic_byrows=.true.)
492 DO WHILE (dbcsr_iterator_blocks_left(iter))
493 CALL dbcsr_iterator_next_block(iter, row, col, data_d, row_offset=rowoff, col_offset=coloff)
494 IF (my_ncol .NE. 1) THEN
495 CALL dgemm('N', 'N', &
496 SIZE(data_d, 1), my_ncol, SIZE(data_d, 2), &
497 my_alpha, data_d(1, 1), SIZE(data_d, 1), &
498 vec_b(coloff, 1), SIZE(vec_b, 1), &
499 1.0_dp, vec_c(rowoff, 1), SIZE(vec_c, 1))
500 ELSE
501 CALL dgemv('N', SIZE(data_d, 1), SIZE(data_d, 2), &
502 my_alpha, data_d(1, 1), SIZE(data_d, 1), &
503 vec_b(coloff, 1), 1, &
504 1.0_dp, vec_c(rowoff, 1), 1)
505 END IF
506 END DO
507 CALL dbcsr_iterator_stop(iter)
508!$OMP END PARALLEL
509
510 ! FIXME ... in the symmetric case, the writes to vec_c depend on the column, not the row. This makes OMP-ing more difficult
511 ! needs e.g. a buffer for vec_c and a reduction of that buffer.
512 IF (has_symm) THEN
513 CALL dbcsr_iterator_readonly_start(iter, matrix_a)
514 DO WHILE (dbcsr_iterator_blocks_left(iter))
515 CALL dbcsr_iterator_next_block(iter, row, col, data_d, row_offset=rowoff, col_offset=coloff)
516 IF (row .NE. col) THEN
517 IF (my_ncol .NE. 1) THEN
518 CALL dgemm('T', 'N', &
519 SIZE(data_d, 2), my_ncol, SIZE(data_d, 1), &
520 my_alpha2, data_d(1, 1), SIZE(data_d, 1), &
521 vec_b(rowoff, 1), SIZE(vec_b, 1), &
522 1.0_dp, vec_c(coloff, 1), SIZE(vec_c, 1))
523 ELSE
524 CALL dgemv('T', SIZE(data_d, 1), SIZE(data_d, 2), &
525 my_alpha2, data_d(1, 1), SIZE(data_d, 1), &
526 vec_b(rowoff, 1), 1, &
527 1.0_dp, vec_c(coloff, 1), 1)
528 END IF
529 END IF
530 END DO
531 CALL dbcsr_iterator_stop(iter)
532 END IF
533
534 CALL timestop(timing_handle)
535 END SUBROUTINE dbcsr_multiply_local
536
537! **************************************************************************************************
538!> \brief multiply a dbcsr with a fm matrix
539!>
540!> For backwards compatibility with BLAS XGEMM, this routine supports
541!> the multiplication of matrices with incompatible dimensions.
542!>
543!> \param[in] matrix DBCSR matrix
544!> \param fm_in full matrix
545!> \param fm_out full matrix
546!> \param[in] ncol nbr of columns
547!> \param[in] alpha alpha
548!> \param[in] beta beta
549!>
550! **************************************************************************************************
551 SUBROUTINE cp_dbcsr_sm_fm_multiply(matrix, fm_in, fm_out, ncol, alpha, beta)
552 TYPE(dbcsr_type), INTENT(IN) :: matrix
553 TYPE(cp_fm_type), INTENT(IN) :: fm_in
554 TYPE(cp_fm_type), INTENT(INOUT) :: fm_out
555 INTEGER, INTENT(IN) :: ncol
556 REAL(dp), INTENT(IN), OPTIONAL :: alpha, beta
557
558 CHARACTER(LEN=*), PARAMETER :: routinen = 'cp_dbcsr_sm_fm_multiply'
559
560 INTEGER :: a_ncol, a_nrow, b_ncol, b_nrow, c_ncol, &
561 c_nrow, k_in, k_out, timing_handle, &
562 timing_handle_mult
563 INTEGER, DIMENSION(:), POINTER :: col_blk_size, col_blk_size_right_in, &
564 col_blk_size_right_out, col_dist, &
565 row_blk_size, row_dist
566 TYPE(dbcsr_type) :: in, out
567 TYPE(dbcsr_distribution_type) :: dist, dist_right_in, product_dist
568 REAL(dp) :: my_alpha, my_beta
569
570 CALL timeset(routinen, timing_handle)
571
572 my_alpha = 1.0_dp
573 my_beta = 0.0_dp
574 IF (PRESENT(alpha)) my_alpha = alpha
575 IF (PRESENT(beta)) my_beta = beta
576
577 ! TODO
578 CALL cp_fm_get_info(fm_in, ncol_global=b_ncol, nrow_global=b_nrow)
579 CALL cp_fm_get_info(fm_out, ncol_global=c_ncol, nrow_global=c_nrow)
580 CALL dbcsr_get_info(matrix, nfullrows_total=a_nrow, nfullcols_total=a_ncol)
581 !WRITE(*,*) "cp_dbcsr_sm_fm_multiply: A ", a_nrow, "x", a_ncol
582 !WRITE(*,*) "cp_dbcsr_sm_fm_multiply: B ", b_nrow, "x", b_ncol
583 !WRITE(*,*) "cp_dbcsr_sm_fm_multiply: C ", c_nrow, "x", c_ncol
584
585 CALL cp_fm_get_info(fm_out, ncol_global=k_out)
586
587 CALL cp_fm_get_info(fm_in, ncol_global=k_in)
588 !write(*,*)routineN//" -----------------------------------"
589 !IF (k_in .NE. k_out) &
590 ! WRITE(*,'(3(A,I5,1X),2(A,F5.2,1X))')&
591 ! routineN//" ncol", ncol,'k_in',k_in,'k_out',k_out,&
592 ! 'alpha',my_alpha,'beta',my_beta
593
594 IF (ncol .GT. 0 .AND. k_out .GT. 0 .AND. k_in .GT. 0) THEN
595 CALL dbcsr_get_info(matrix, row_blk_size=row_blk_size, col_blk_size=col_blk_size, distribution=dist)
596 CALL dbcsr_create_dist_r_unrot(dist_right_in, dist, k_in, col_blk_size_right_in)
597
598 CALL dbcsr_create(in, "D", dist_right_in, dbcsr_type_no_symmetry, &
599 col_blk_size, col_blk_size_right_in)
600
601 CALL dbcsr_distribution_get(dist, row_dist=row_dist)
602 CALL dbcsr_distribution_get(dist_right_in, col_dist=col_dist)
603 CALL dbcsr_distribution_new(product_dist, template=dist, &
604 row_dist=row_dist, col_dist=col_dist)
605 ALLOCATE (col_blk_size_right_out(SIZE(col_blk_size_right_in)))
606 col_blk_size_right_out = col_blk_size_right_in
607 CALL match_col_sizes(col_blk_size_right_out, col_blk_size_right_in, k_out)
608
609 !if (k_in .ne. k_out) then
610 ! write(*,*)routineN//" in cs", col_blk_size_right_in
611 ! write(*,*)routineN//" out cs", col_blk_size_right_out
612 !endif
613
614 CALL dbcsr_create(out, "D", product_dist, dbcsr_type_no_symmetry, &
615 row_blk_size, col_blk_size_right_out)
616
617 CALL copy_fm_to_dbcsr(fm_in, in)
618 IF (ncol .NE. k_out .OR. my_beta .NE. 0.0_dp) &
619 CALL copy_fm_to_dbcsr(fm_out, out)
620
621 CALL timeset(routinen//'_core', timing_handle_mult)
622 CALL dbcsr_multiply("N", "N", my_alpha, matrix, in, my_beta, out, &
623 last_column=ncol)
624 CALL timestop(timing_handle_mult)
625
626 CALL copy_dbcsr_to_fm(out, fm_out)
627
628 CALL dbcsr_release(in)
629 CALL dbcsr_release(out)
630 DEALLOCATE (col_blk_size_right_in, col_blk_size_right_out)
631 CALL dbcsr_distribution_release(dist_right_in)
632 CALL dbcsr_distribution_release(product_dist)
633
634 END IF
635
636 CALL timestop(timing_handle)
637
638 END SUBROUTINE cp_dbcsr_sm_fm_multiply
639
640! **************************************************************************************************
641!> \brief ...
642!> \param sizes1 ...
643!> \param sizes2 ...
644!> \param full_num ...
645! **************************************************************************************************
646 SUBROUTINE match_col_sizes(sizes1, sizes2, full_num)
647 INTEGER, DIMENSION(:), INTENT(INOUT) :: sizes1
648 INTEGER, DIMENSION(:), INTENT(IN) :: sizes2
649 INTEGER, INTENT(IN) :: full_num
650
651 INTEGER :: left, n1, n2, p, rm, used
652
653 n1 = SIZE(sizes1)
654 n2 = SIZE(sizes2)
655 IF (n1 .NE. n2) &
656 cpabort("distributions must be equal!")
657 sizes1(1:n1) = sizes2(1:n1)
658 used = sum(sizes1(1:n1))
659 ! If sizes1 does not cover everything, then we increase the
660 ! size of the last block; otherwise we reduce the blocks
661 ! (from the end) until it is small enough.
662 IF (used .LT. full_num) THEN
663 sizes1(n1) = sizes1(n1) + full_num - used
664 ELSE
665 left = used - full_num
666 p = n1
667 DO WHILE (left .GT. 0 .AND. p .GT. 0)
668 rm = min(left, sizes1(p))
669 sizes1(p) = sizes1(p) - rm
670 left = left - rm
671 p = p - 1
672 END DO
673 END IF
674 END SUBROUTINE match_col_sizes
675
676! **************************************************************************************************
677!> \brief performs the multiplication sparse_matrix+dense_mat*dens_mat^T
678!> if matrix_g is not explicitly given, matrix_v^T will be used
679!> this can be important to save the necessary redistribute for a
680!> different matrix_g and increase performance.
681!> \param sparse_matrix ...
682!> \param matrix_v ...
683!> \param matrix_g ...
684!> \param ncol ...
685!> \param alpha ...
686!> \param keep_sparsity Determines if the sparsity of sparse_matrix is retained
687!> by default it is TRUE
688!> \param symmetry_mode There are the following modes
689!> 1: sparse_matrix += 0.5*alpha*(v*g^T+g^T*v) (symmetric update)
690!> -1: sparse_matrix += 0.5*alpha*(v*g^T-g^T*v) (skewsymmetric update)
691!> else: sparse_matrix += alpha*v*g^T (no symmetry, default)
692!> saves some redistribution steps
693! **************************************************************************************************
694 SUBROUTINE cp_dbcsr_plus_fm_fm_t(sparse_matrix, matrix_v, matrix_g, ncol, alpha, keep_sparsity, symmetry_mode)
695 TYPE(dbcsr_type), INTENT(INOUT) :: sparse_matrix
696 TYPE(cp_fm_type), INTENT(IN) :: matrix_v
697 TYPE(cp_fm_type), INTENT(IN), OPTIONAL :: matrix_g
698 INTEGER, INTENT(IN) :: ncol
699 REAL(kind=dp), INTENT(IN), OPTIONAL :: alpha
700 LOGICAL, INTENT(IN), OPTIONAL :: keep_sparsity
701 INTEGER, INTENT(IN), OPTIONAL :: symmetry_mode
702
703 CHARACTER(LEN=*), PARAMETER :: routinen = 'cp_dbcsr_plus_fm_fm_t'
704
705 INTEGER :: k, my_symmetry_mode, nao, npcols, &
706 timing_handle
707 INTEGER, DIMENSION(:), POINTER :: col_blk_size_left, col_dist_left, &
708 row_blk_size, row_dist
709 LOGICAL :: check_product, my_keep_sparsity
710 REAL(kind=dp) :: my_alpha, norm
711 TYPE(cp_fm_struct_type), POINTER :: fm_struct_tmp
712 TYPE(cp_fm_type) :: fm_matrix
713 TYPE(dbcsr_distribution_type) :: dist_left, sparse_dist
714 TYPE(dbcsr_type) :: mat_g, mat_v, sparse_matrix2, &
715 sparse_matrix3
716
717 check_product = .false.
718
719 CALL timeset(routinen, timing_handle)
720
721 my_keep_sparsity = .true.
722 IF (PRESENT(keep_sparsity)) my_keep_sparsity = keep_sparsity
723
724 my_symmetry_mode = 0
725 IF (PRESENT(symmetry_mode)) my_symmetry_mode = symmetry_mode
726
727 NULLIFY (col_dist_left)
728
729 IF (ncol .GT. 0) THEN
730 IF (.NOT. dbcsr_valid_index(sparse_matrix)) &
731 cpabort("sparse_matrix must pre-exist")
732 !
733 ! Setup matrix_v
734 CALL cp_fm_get_info(matrix_v, ncol_global=k)
735 !WRITE(*,*)routineN//'truncated mult k, ncol',k,ncol,' PRESENT (matrix_g)',PRESENT (matrix_g)
736 CALL dbcsr_get_info(sparse_matrix, distribution=sparse_dist)
737 CALL dbcsr_distribution_get(sparse_dist, npcols=npcols, row_dist=row_dist)
738 CALL create_bl_distribution(col_dist_left, col_blk_size_left, k, npcols)
739 CALL dbcsr_distribution_new(dist_left, template=sparse_dist, &
740 row_dist=row_dist, col_dist=col_dist_left)
741 DEALLOCATE (col_dist_left)
742 CALL dbcsr_get_info(sparse_matrix, row_blk_size=row_blk_size)
743 CALL dbcsr_create(mat_v, "DBCSR matrix_v", dist_left, dbcsr_type_no_symmetry, &
744 row_blk_size, col_blk_size_left)
745 CALL copy_fm_to_dbcsr(matrix_v, mat_v)
746 CALL dbcsr_verify_matrix(mat_v)
747 !
748 ! Setup matrix_g
749 IF (PRESENT(matrix_g)) THEN
750 CALL dbcsr_create(mat_g, "DBCSR matrix_g", dist_left, dbcsr_type_no_symmetry, &
751 row_blk_size, col_blk_size_left)
752 CALL copy_fm_to_dbcsr(matrix_g, mat_g)
753 END IF
754 !
755 DEALLOCATE (col_blk_size_left)
756 CALL dbcsr_distribution_release(dist_left)
757 !
758 !
759 IF (check_product) THEN
760 CALL cp_fm_get_info(matrix_v, nrow_global=nao)
761 CALL cp_fm_struct_create(fm_struct_tmp, context=matrix_v%matrix_struct%context, nrow_global=nao, &
762 ncol_global=nao, para_env=matrix_v%matrix_struct%para_env)
763 CALL cp_fm_create(fm_matrix, fm_struct_tmp, name="fm matrix")
764 CALL cp_fm_struct_release(fm_struct_tmp)
765 CALL copy_dbcsr_to_fm(sparse_matrix, fm_matrix)
766 CALL dbcsr_copy(sparse_matrix3, sparse_matrix)
767 END IF
768 !
769 my_alpha = 1.0_dp
770 IF (PRESENT(alpha)) my_alpha = alpha
771 IF (PRESENT(matrix_g)) THEN
772 IF (my_symmetry_mode == 1) THEN
773 ! Symmetric mode
774 CALL dbcsr_multiply("N", "T", 0.5_dp*my_alpha, mat_v, mat_g, &
775 1.0_dp, sparse_matrix, &
776 retain_sparsity=my_keep_sparsity, &
777 last_k=ncol)
778 CALL dbcsr_multiply("N", "T", 0.5_dp*my_alpha, mat_g, mat_v, &
779 1.0_dp, sparse_matrix, &
780 retain_sparsity=my_keep_sparsity, &
781 last_k=ncol)
782 ELSE IF (my_symmetry_mode == -1) THEN
783 ! Skewsymmetric mode
784 CALL dbcsr_multiply("N", "T", 0.5_dp*my_alpha, mat_v, mat_g, &
785 1.0_dp, sparse_matrix, &
786 retain_sparsity=my_keep_sparsity, &
787 last_k=ncol)
788 CALL dbcsr_multiply("N", "T", -0.5_dp*my_alpha, mat_g, mat_v, &
789 1.0_dp, sparse_matrix, &
790 retain_sparsity=my_keep_sparsity, &
791 last_k=ncol)
792 ELSE
793 ! Normal mode
794 CALL dbcsr_multiply("N", "T", my_alpha, mat_v, mat_g, &
795 1.0_dp, sparse_matrix, &
796 retain_sparsity=my_keep_sparsity, &
797 last_k=ncol)
798 END IF
799 ELSE
800 CALL dbcsr_multiply("N", "T", my_alpha, mat_v, mat_v, &
801 1.0_dp, sparse_matrix, &
802 retain_sparsity=my_keep_sparsity, &
803 last_k=ncol)
804 END IF
805
806 IF (check_product) THEN
807 IF (PRESENT(matrix_g)) THEN
808 IF (my_symmetry_mode == 1) THEN
809 CALL cp_fm_gemm("N", "T", nao, nao, ncol, 0.5_dp*my_alpha, matrix_v, matrix_g, &
810 1.0_dp, fm_matrix)
811 CALL cp_fm_gemm("N", "T", nao, nao, ncol, 0.5_dp*my_alpha, matrix_g, matrix_v, &
812 1.0_dp, fm_matrix)
813 ELSE IF (my_symmetry_mode == -1) THEN
814 CALL cp_fm_gemm("N", "T", nao, nao, ncol, 0.5_dp*my_alpha, matrix_v, matrix_g, &
815 1.0_dp, fm_matrix)
816 CALL cp_fm_gemm("N", "T", nao, nao, ncol, -0.5_dp*my_alpha, matrix_g, matrix_v, &
817 1.0_dp, fm_matrix)
818 ELSE
819 CALL cp_fm_gemm("N", "T", nao, nao, ncol, my_alpha, matrix_v, matrix_g, &
820 1.0_dp, fm_matrix)
821 END IF
822 ELSE
823 CALL cp_fm_gemm("N", "T", nao, nao, ncol, my_alpha, matrix_v, matrix_v, &
824 1.0_dp, fm_matrix)
825 END IF
826
827 CALL dbcsr_copy(sparse_matrix2, sparse_matrix)
828 CALL dbcsr_scale(sparse_matrix2, alpha_scalar=0.0_dp)
829 CALL copy_fm_to_dbcsr(fm_matrix, sparse_matrix2, keep_sparsity=my_keep_sparsity)
830 CALL dbcsr_add(sparse_matrix2, sparse_matrix, alpha_scalar=1.0_dp, &
831 beta_scalar=-1.0_dp)
832 norm = dbcsr_frobenius_norm(sparse_matrix2)
833 WRITE (*, *) 'nao=', nao, ' k=', k, ' ncol=', ncol, ' my_alpha=', my_alpha
834 WRITE (*, *) 'PRESENT (matrix_g)', PRESENT(matrix_g)
835 WRITE (*, *) 'matrix_type=', dbcsr_get_matrix_type(sparse_matrix)
836 WRITE (*, *) 'norm(sm+alpha*v*g^t - fm+alpha*v*g^t)/n=', norm/real(nao, dp)
837 IF (norm/real(nao, dp) .GT. 1e-12_dp) THEN
838 !WRITE(*,*) 'fm_matrix'
839 !DO j=1,SIZE(fm_matrix%local_data,2)
840 ! DO i=1,SIZE(fm_matrix%local_data,1)
841 ! WRITE(*,'(A,I3,A,I3,A,E26.16,A)') 'a(',i,',',j,')=',fm_matrix%local_data(i,j),';'
842 ! ENDDO
843 !ENDDO
844 !WRITE(*,*) 'mat_v'
845 !CALL dbcsr_print(mat_v)
846 !WRITE(*,*) 'mat_g'
847 !CALL dbcsr_print(mat_g)
848 !WRITE(*,*) 'sparse_matrix'
849 !CALL dbcsr_print(sparse_matrix)
850 !WRITE(*,*) 'sparse_matrix2 (-sm + sparse(fm))'
851 !CALL dbcsr_print(sparse_matrix2)
852 !WRITE(*,*) 'sparse_matrix3 (copy of sm input)'
853 !CALL dbcsr_print(sparse_matrix3)
854 !stop
855 END IF
856 CALL dbcsr_release(sparse_matrix2)
857 CALL dbcsr_release(sparse_matrix3)
858 CALL cp_fm_release(fm_matrix)
859 END IF
860 CALL dbcsr_release(mat_v)
861 IF (PRESENT(matrix_g)) CALL dbcsr_release(mat_g)
862 END IF
863 CALL timestop(timing_handle)
864
865 END SUBROUTINE cp_dbcsr_plus_fm_fm_t
866
867! **************************************************************************************************
868!> \brief Utility function to copy a specially shaped fm to dbcsr_matrix
869!> The result matrix will be the matrix in dbcsr format
870!> with the row blocks sizes according to the block_sizes of the template
871!> and the col blocks sizes evenly blocked with the internal dbcsr conversion
872!> size (32 is the current default)
873!> \param matrix ...
874!> \param fm_in ...
875!> \param template ...
876! **************************************************************************************************
877 SUBROUTINE cp_fm_to_dbcsr_row_template(matrix, fm_in, template)
878 TYPE(dbcsr_type), INTENT(INOUT) :: matrix
879 TYPE(cp_fm_type), INTENT(IN) :: fm_in
880 TYPE(dbcsr_type), INTENT(IN) :: template
881
882 INTEGER :: k_in
883 INTEGER, DIMENSION(:), POINTER :: col_blk_size_right_in, row_blk_size
884 TYPE(dbcsr_distribution_type) :: dist_right_in, tmpl_dist
885
886 CALL cp_fm_get_info(fm_in, ncol_global=k_in)
887
888 CALL dbcsr_get_info(template, distribution=tmpl_dist)
889 CALL dbcsr_create_dist_r_unrot(dist_right_in, tmpl_dist, k_in, col_blk_size_right_in)
890 CALL dbcsr_get_info(template, row_blk_size=row_blk_size)
891 CALL dbcsr_create(matrix, "D", dist_right_in, dbcsr_type_no_symmetry, &
892 row_blk_size, col_blk_size_right_in)
893
894 CALL copy_fm_to_dbcsr(fm_in, matrix)
895 DEALLOCATE (col_blk_size_right_in)
896 CALL dbcsr_distribution_release(dist_right_in)
897
898 END SUBROUTINE cp_fm_to_dbcsr_row_template
899
900! **************************************************************************************************
901!> \brief Utility function to create an arbitrary shaped dbcsr matrix
902!> with the same processor grid as the template matrix
903!> both row sizes and col sizes are evenly blocked with the internal
904!> dbcsr_conversion size (32 is the current default)
905!> \param matrix dbcsr matrix to be created
906!> \param template template dbcsr matrix giving its mp_env
907!> \param m global row size of output matrix
908!> \param n global col size of output matrix
909!> \param sym ...
910! **************************************************************************************************
911 SUBROUTINE cp_dbcsr_m_by_n_from_template(matrix, template, m, n, sym)
912 TYPE(dbcsr_type), INTENT(INOUT) :: matrix, template
913 INTEGER, INTENT(IN) :: m, n
914 CHARACTER, INTENT(IN), OPTIONAL :: sym
915
916 CHARACTER :: mysym
917 INTEGER :: npcols, nprows
918 INTEGER, DIMENSION(:), POINTER :: col_blk_size, col_dist, row_blk_size, &
919 row_dist
920 TYPE(dbcsr_distribution_type) :: dist_m_n, tmpl_dist
921
922 CALL dbcsr_get_info(template, matrix_type=mysym, distribution=tmpl_dist)
923
924 IF (PRESENT(sym)) mysym = sym
925
926 NULLIFY (row_dist, col_dist)
927 NULLIFY (row_blk_size, col_blk_size)
928 !NULLIFY (row_cluster, col_cluster)
929
930 CALL dbcsr_distribution_get(tmpl_dist, nprows=nprows, npcols=npcols)
931 CALL create_bl_distribution(row_dist, row_blk_size, m, nprows)
932 CALL create_bl_distribution(col_dist, col_blk_size, n, npcols)
933 CALL dbcsr_distribution_new(dist_m_n, template=tmpl_dist, &
934 row_dist=row_dist, col_dist=col_dist, &
935 !row_cluster=row_cluster, col_cluster=col_cluster, &
936 reuse_arrays=.true.)
937
938 CALL dbcsr_create(matrix, "m_n_template", dist_m_n, mysym, &
939 row_blk_size, col_blk_size, reuse_arrays=.true.)
940 CALL dbcsr_distribution_release(dist_m_n)
941
942 END SUBROUTINE cp_dbcsr_m_by_n_from_template
943
944! **************************************************************************************************
945!> \brief Utility function to create dbcsr matrix, m x n matrix (n arbitrary)
946!> with the same processor grid and row distribution as the template matrix
947!> col sizes are evenly blocked with the internal
948!> dbcsr_conversion size (32 is the current default)
949!> \param matrix dbcsr matrix to be created
950!> \param template template dbcsr matrix giving its mp_env
951!> \param n global col size of output matrix
952!> \param sym ...
953! **************************************************************************************************
954 SUBROUTINE cp_dbcsr_m_by_n_from_row_template(matrix, template, n, sym)
955 TYPE(dbcsr_type), INTENT(INOUT) :: matrix, template
956 INTEGER :: n
957 CHARACTER, OPTIONAL :: sym
958
959 CHARACTER :: mysym
960 INTEGER :: npcols
961 INTEGER, DIMENSION(:), POINTER :: col_blk_size, col_dist, row_blk_size, &
962 row_dist
963 TYPE(dbcsr_distribution_type) :: dist_m_n, tmpl_dist
964
965 mysym = dbcsr_get_matrix_type(template)
966 IF (PRESENT(sym)) mysym = sym
967
968 CALL dbcsr_get_info(template, distribution=tmpl_dist)
969 CALL dbcsr_distribution_get(tmpl_dist, &
970 npcols=npcols, &
971 row_dist=row_dist)
972
973 NULLIFY (col_dist, col_blk_size)
974 CALL create_bl_distribution(col_dist, col_blk_size, n, npcols)
975 CALL dbcsr_distribution_new(dist_m_n, template=tmpl_dist, &
976 row_dist=row_dist, col_dist=col_dist)
977
978 CALL dbcsr_get_info(template, row_blk_size=row_blk_size)
979 CALL dbcsr_create(matrix, "m_n_template", dist_m_n, mysym, row_blk_size, col_blk_size)
980
981 DEALLOCATE (col_dist, col_blk_size)
982 CALL dbcsr_distribution_release(dist_m_n)
983
984 END SUBROUTINE cp_dbcsr_m_by_n_from_row_template
985
986! **************************************************************************************************
987!> \brief Distributes elements into blocks and into bins
988!>
989!> \param[out] block_distribution block distribution to bins
990!> \param[out] block_size sizes of blocks
991!> \param[in] nelements number of elements to bin
992!> \param[in] nbins number of bins
993!> \par Term clarification
994!> An example: blocks are atom blocks and bins are process rows/columns.
995! **************************************************************************************************
996 SUBROUTINE create_bl_distribution(block_distribution, &
997 block_size, nelements, nbins)
998 INTEGER, DIMENSION(:), INTENT(OUT), POINTER :: block_distribution, block_size
999 INTEGER, INTENT(IN) :: nelements, nbins
1000
1001 CHARACTER(len=*), PARAMETER :: routinen = 'create_bl_distribution', &
1002 routinep = modulen//':'//routinen
1003
1004 INTEGER :: bin, blk_layer, element_stack, els, &
1005 estimated_blocks, max_blocks_per_bin, &
1006 nblks, nblocks, stat
1007 INTEGER, DIMENSION(:), POINTER :: blk_dist, blk_sizes
1008
1009! ---------------------------------------------------------------------------
1010
1011 NULLIFY (block_distribution)
1012 NULLIFY (block_size)
1013 ! Define the sizes on which we build the distribution.
1014 IF (nelements .GT. 0) THEN
1015
1016 nblocks = ceiling(real(nelements, kind=dp)/real(max_elements_per_block, kind=dp))
1017 max_blocks_per_bin = ceiling(real(nblocks, kind=dp)/real(nbins, kind=dp))
1018
1019 IF (debug_mod) THEN
1020 WRITE (*, '(1X,A,1X,A,I7,A,I7,A)') routinep, "For", nelements, &
1021 " elements and", nbins, " bins"
1022 WRITE (*, '(1X,A,1X,A,I7,A)') routinep, "There are", &
1023 max_elements_per_block, " max elements per block"
1024 WRITE (*, '(1X,A,1X,A,I7,A)') routinep, "There are", &
1025 nblocks, " blocks"
1026 WRITE (*, '(1X,A,1X,A,I7,A)') routinep, "There are", &
1027 max_blocks_per_bin, " max blocks/bin"
1028 END IF
1029
1030 estimated_blocks = max_blocks_per_bin*nbins
1031 ALLOCATE (blk_dist(estimated_blocks), stat=stat)
1032 IF (stat /= 0) &
1033 cpabort("blk_dist")
1034 ALLOCATE (blk_sizes(estimated_blocks), stat=stat)
1035 IF (stat /= 0) &
1036 cpabort("blk_sizes")
1037 element_stack = 0
1038 nblks = 0
1039 DO blk_layer = 1, max_blocks_per_bin
1040 DO bin = 0, nbins - 1
1041 els = min(max_elements_per_block, nelements - element_stack)
1042 IF (els .GT. 0) THEN
1043 element_stack = element_stack + els
1044 nblks = nblks + 1
1045 blk_dist(nblks) = bin
1046 blk_sizes(nblks) = els
1047 IF (debug_mod) WRITE (*, '(1X,A,I5,A,I5,A,I5)') routinep//" Assigning", &
1048 els, " elements as block", nblks, " to bin", bin
1049 END IF
1050 END DO
1051 END DO
1052 ! Create the output arrays.
1053 IF (nblks .EQ. estimated_blocks) THEN
1054 block_distribution => blk_dist
1055 block_size => blk_sizes
1056 ELSE
1057 ALLOCATE (block_distribution(nblks), stat=stat)
1058 IF (stat /= 0) &
1059 cpabort("blk_dist")
1060 block_distribution(:) = blk_dist(1:nblks)
1061 DEALLOCATE (blk_dist)
1062 ALLOCATE (block_size(nblks), stat=stat)
1063 IF (stat /= 0) &
1064 cpabort("blk_sizes")
1065 block_size(:) = blk_sizes(1:nblks)
1066 DEALLOCATE (blk_sizes)
1067 END IF
1068 ELSE
1069 ALLOCATE (block_distribution(0), stat=stat)
1070 IF (stat /= 0) &
1071 cpabort("blk_dist")
1072 ALLOCATE (block_size(0), stat=stat)
1073 IF (stat /= 0) &
1074 cpabort("blk_sizes")
1075 END IF
10761579 FORMAT(i5, 1x, i5, 1x, i5, 1x, i5, 1x, i5, 1x, i5, 1x, i5, 1x, i5, 1x, i5, 1x, i5)
1077 IF (debug_mod) THEN
1078 WRITE (*, '(1X,A,A)') routinep//" Distribution"
1079 WRITE (*, 1579) block_distribution(:)
1080 WRITE (*, '(1X,A,A)') routinep//" Sizes"
1081 WRITE (*, 1579) block_size(:)
1082 END IF
1083 END SUBROUTINE create_bl_distribution
1084
1085! **************************************************************************************************
1086!> \brief Creates a new distribution for the right matrix in a matrix
1087!> multiplication with unrotated grid.
1088!> \param[out] dist_right new distribution for the right matrix
1089!> \param[in] dist_left the distribution of the left matrix
1090!> \param[in] ncolumns number of columns in right matrix
1091!> \param[out] right_col_blk_sizes sizes of blocks in the created column
1092!> \par The new row distribution for the right matrix is the same as the row
1093!> distribution of the left matrix, while the column distribution is
1094!> created so that it is appropriate to the parallel environment.
1095! **************************************************************************************************
1096 SUBROUTINE dbcsr_create_dist_r_unrot(dist_right, dist_left, ncolumns, &
1097 right_col_blk_sizes)
1098 TYPE(dbcsr_distribution_type), INTENT(OUT) :: dist_right
1099 TYPE(dbcsr_distribution_type), INTENT(IN) :: dist_left
1100 INTEGER, INTENT(IN) :: ncolumns
1101 INTEGER, DIMENSION(:), INTENT(OUT), POINTER :: right_col_blk_sizes
1102
1103 INTEGER :: multiplicity, ncols, nimages, npcols, &
1104 nprows
1105 INTEGER, ALLOCATABLE, DIMENSION(:) :: tmp_images
1106 INTEGER, DIMENSION(:), POINTER :: old_col_dist, right_col_dist, &
1107 right_row_dist
1108
1109 CALL dbcsr_distribution_get(dist_left, &
1110 ncols=ncols, &
1111 col_dist=old_col_dist, &
1112 nprows=nprows, &
1113 npcols=npcols)
1114
1115 ! Create the column distribution
1116 CALL create_bl_distribution(right_col_dist, right_col_blk_sizes, ncolumns, npcols)
1117 ! Create an even row distribution.
1118 ALLOCATE (right_row_dist(ncols), tmp_images(ncols))
1119 nimages = lcm(nprows, npcols)/nprows
1120 multiplicity = nprows/gcd(nprows, npcols)
1121 CALL rebin_distribution(right_row_dist, tmp_images, old_col_dist, nprows, multiplicity, nimages)
1122
1123 CALL dbcsr_distribution_new(dist_right, &
1124 template=dist_left, &
1125 row_dist=right_row_dist, &
1126 col_dist=right_col_dist, &
1127 !row_cluster=dummy,&
1128 !col_cluster=dummy,&
1129 reuse_arrays=.true.)
1130 DEALLOCATE (tmp_images)
1131 END SUBROUTINE dbcsr_create_dist_r_unrot
1132
1133! **************************************************************************************************
1134!> \brief Makes new distribution with decimation and multiplicity
1135!> \param[out] new_bins new real distribution
1136!> \param[out] images new image distribution
1137!> \param[in] source_bins Basis for the new distribution and images
1138!> \param[in] nbins number of bins in the new real distribution
1139!> \param[in] multiplicity multiplicity
1140!> \param[in] nimages number of images in the new distribution
1141!> \par Definition of multiplicity and nimages
1142!> Multiplicity and decimation (number of images) are used to
1143!> match process grid coordinates on non-square process
1144!> grids. Given source_nbins and target_nbins, their relation is
1145!> source_nbins * target_multiplicity
1146!> = target_nbins * target_nimages.
1147!> It is best when both multiplicity and nimages are small. To
1148!> get these two factors, then, one can use the following formulas:
1149!> nimages = lcm(source_nbins, target_nbins) / target_nbins
1150!> multiplicity = target_nbins / gcd(source_nbins, target_nbins)
1151!> from the target's point of view (nimages = target_nimages).
1152!> \par Mapping
1153!> The new distribution comprises of real bins and images within
1154!> bins. These can be view as target_nbins*nimages virtual
1155!> columns. These same virtual columns are also
1156!> source_nbins*multiplicity in number. Therefore these virtual
1157!> columns are mapped from source_nbins*multiplicity onto
1158!> target_bins*nimages (each target bin has nimages images):
1159!> Source 4: |1 2 3|4 5 6|7 8 9|A B C| (4*3)
1160!> Target 6: |1 2|3 4|5 6|7 8|9 A|B C| (6*2)
1161!> multiplicity=3, nimages=2, 12 virtual columns (1-C).
1162!> Source bin elements are evenly mapped into one of multiplicity
1163!> virtual columns. Other (non-even, block-size aware) mappings
1164!> could be better.
1165! **************************************************************************************************
1166 SUBROUTINE rebin_distribution(new_bins, images, source_bins, &
1167 nbins, multiplicity, nimages)
1168 INTEGER, DIMENSION(:), INTENT(OUT) :: new_bins, images
1169 INTEGER, DIMENSION(:), INTENT(IN) :: source_bins
1170 INTEGER, INTENT(IN) :: nbins, multiplicity, nimages
1171
1172 INTEGER :: bin, i, old_nbins, virtual_bin
1173 INTEGER, ALLOCATABLE, DIMENSION(:) :: bin_multiplier
1174
1175! ---------------------------------------------------------------------------
1176
1177 IF (mod(nbins*nimages, multiplicity) .NE. 0) THEN
1178 cpwarn("mulitplicity is not divisor of new process grid coordinate")
1179 END IF
1180 old_nbins = (nbins*nimages)/multiplicity
1181 ALLOCATE (bin_multiplier(0:old_nbins - 1))
1182 bin_multiplier(:) = 0
1183 DO i = 1, SIZE(new_bins)
1184 IF (i .LE. SIZE(source_bins)) THEN
1185 bin = source_bins(i)
1186 ELSE
1187 ! Fill remainder with a cyclic distribution
1188 bin = mod(i, old_nbins)
1189 END IF
1190 virtual_bin = bin*multiplicity + bin_multiplier(bin)
1191 new_bins(i) = virtual_bin/nimages
1192 images(i) = 1 + mod(virtual_bin, nimages)
1193 bin_multiplier(bin) = bin_multiplier(bin) + 1
1194 IF (bin_multiplier(bin) .GE. multiplicity) THEN
1195 bin_multiplier(bin) = 0
1196 END IF
1197 END DO
1198 END SUBROUTINE rebin_distribution
1199
1200! **************************************************************************************************
1201!> \brief Creates a block-cyclic compatible distribution
1202!>
1203!> All blocks in a dimension, except for possibly the last
1204!> block, have the same size.
1205!> \param[out] dist the elemental distribution
1206!> \param[in] nrows number of full rows
1207!> \param[in] ncolumns number of full columns
1208!> \param[in] nrow_block size of row blocks
1209!> \param[in] ncol_block size of column blocks
1210!> \param group_handle ...
1211!> \param pgrid ...
1212!> \param[out] row_blk_sizes row block sizes
1213!> \param[out] col_blk_sizes column block sizes
1214! **************************************************************************************************
1215 SUBROUTINE dbcsr_create_dist_block_cyclic(dist, nrows, ncolumns, &
1216 nrow_block, ncol_block, group_handle, pgrid, row_blk_sizes, col_blk_sizes)
1217 TYPE(dbcsr_distribution_type), INTENT(OUT) :: dist
1218 INTEGER, INTENT(IN) :: nrows, ncolumns, nrow_block, ncol_block, &
1219 group_handle
1220 INTEGER, DIMENSION(:, :), POINTER :: pgrid
1221 INTEGER, DIMENSION(:), INTENT(OUT), POINTER :: row_blk_sizes, col_blk_sizes
1222
1223 CHARACTER(len=*), PARAMETER :: routinen = 'dbcsr_create_dist_block_cyclic'
1224
1225 INTEGER :: nblkcols, nblkrows, npcols, nprows, &
1226 pdim, sz
1227 INTEGER, DIMENSION(:), POINTER :: cd_data, rd_data
1228
1229 ! Row sizes
1230 IF (nrow_block .EQ. 0) THEN
1231 nblkrows = 0
1232 sz = 0
1233 ELSE
1234 nblkrows = nrows/nrow_block
1235 sz = mod(nrows, nrow_block)
1236 END IF
1237 IF (sz .GT. 0) nblkrows = nblkrows + 1
1238 ALLOCATE (row_blk_sizes(nblkrows), rd_data(nblkrows))
1239 row_blk_sizes = nrow_block
1240 IF (sz .NE. 0) row_blk_sizes(nblkrows) = sz
1241
1242 ! Column sizes
1243 IF (ncol_block .EQ. 0) THEN
1244 nblkcols = 0
1245 sz = 0
1246 ELSE
1247 nblkcols = ncolumns/ncol_block
1248 sz = mod(ncolumns, ncol_block)
1249 END IF
1250 IF (sz .GT. 0) nblkcols = nblkcols + 1
1251 ALLOCATE (col_blk_sizes(nblkcols), cd_data(nblkcols))
1252 col_blk_sizes = ncol_block
1253 IF (sz .NE. 0) col_blk_sizes(nblkcols) = sz
1254 !
1255 IF (debug_mod) THEN
1256 WRITE (*, *) routinen//" nrows,nrow_block,nblkrows=", &
1257 nrows, nrow_block, nblkrows
1258 WRITE (*, *) routinen//" ncols,ncol_block,nblkcols=", &
1259 ncolumns, ncol_block, nblkcols
1260 END IF
1261 ! Calculate process row distribution
1262 nprows = SIZE(pgrid, 1)
1263 DO pdim = 0, min(nprows - 1, nblkrows - 1)
1264 rd_data(1 + pdim:nblkrows:nprows) = pdim
1265 END DO
1266 ! Calculate process column distribution
1267 npcols = SIZE(pgrid, 2)
1268 DO pdim = 0, min(npcols - 1, nblkcols - 1)
1269 cd_data(1 + pdim:nblkcols:npcols) = pdim
1270 END DO
1271 !
1272 IF (debug_mod) THEN
1273 WRITE (*, *) routinen//" row_dist", &
1274 rd_data
1275 WRITE (*, *) routinen//" col_dist", &
1276 cd_data
1277 END IF
1278 !
1279 CALL dbcsr_distribution_new(dist, &
1280 group=group_handle, pgrid=pgrid, &
1281 row_dist=rd_data, &
1282 col_dist=cd_data, &
1283 reuse_arrays=.true.)
1284
1285 END SUBROUTINE dbcsr_create_dist_block_cyclic
1286
1287! **************************************************************************************************
1288!> \brief Allocate and initialize a real matrix 1-dimensional set.
1289!> \param[in,out] matrix_set Set containing the DBCSR matrices
1290!> \param[in] nmatrix Size of set
1291!> \par History
1292!> 2009-08-17 Adapted from sparse_matrix_type for DBCSR
1293! **************************************************************************************************
1294 SUBROUTINE allocate_dbcsr_matrix_set_1d(matrix_set, nmatrix)
1295 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_set
1296 INTEGER, INTENT(IN) :: nmatrix
1297
1298 INTEGER :: imatrix
1299
1300 IF (ASSOCIATED(matrix_set)) CALL dbcsr_deallocate_matrix_set(matrix_set)
1301 ALLOCATE (matrix_set(nmatrix))
1302 DO imatrix = 1, nmatrix
1303 NULLIFY (matrix_set(imatrix)%matrix)
1304 END DO
1305 END SUBROUTINE allocate_dbcsr_matrix_set_1d
1306
1307! **************************************************************************************************
1308!> \brief Allocate and initialize a real matrix 2-dimensional set.
1309!> \param[in,out] matrix_set Set containing the DBCSR matrix pointer type
1310!> \param[in] nmatrix Size of set
1311!> \param mmatrix ...
1312!> \par History
1313!> 2009-08-17 Adapted from sparse_matrix_type for DBCSR
1314! **************************************************************************************************
1315 SUBROUTINE allocate_dbcsr_matrix_set_2d(matrix_set, nmatrix, mmatrix)
1316 TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: matrix_set
1317 INTEGER, INTENT(IN) :: nmatrix, mmatrix
1318
1319 INTEGER :: imatrix, jmatrix
1320
1321 IF (ASSOCIATED(matrix_set)) CALL dbcsr_deallocate_matrix_set(matrix_set)
1322 ALLOCATE (matrix_set(nmatrix, mmatrix))
1323 DO jmatrix = 1, mmatrix
1324 DO imatrix = 1, nmatrix
1325 NULLIFY (matrix_set(imatrix, jmatrix)%matrix)
1326 END DO
1327 END DO
1328 END SUBROUTINE allocate_dbcsr_matrix_set_2d
1329
1330! **************************************************************************************************
1331!> \brief Allocate and initialize a real matrix 3-dimensional set.
1332!> \param[in,out] matrix_set Set containing the DBCSR matrix pointer type
1333!> \param[in] nmatrix Size of set
1334!> \param mmatrix ...
1335!> \param pmatrix ...
1336!> \par History
1337!> 2009-08-17 Adapted from sparse_matrix_type for DBCSR
1338! **************************************************************************************************
1339 SUBROUTINE allocate_dbcsr_matrix_set_3d(matrix_set, nmatrix, mmatrix, pmatrix)
1340 TYPE(dbcsr_p_type), DIMENSION(:, :, :), POINTER :: matrix_set
1341 INTEGER, INTENT(IN) :: nmatrix, mmatrix, pmatrix
1342
1343 INTEGER :: imatrix, jmatrix, kmatrix
1344
1345 IF (ASSOCIATED(matrix_set)) CALL dbcsr_deallocate_matrix_set(matrix_set)
1346 ALLOCATE (matrix_set(nmatrix, mmatrix, pmatrix))
1347 DO kmatrix = 1, pmatrix
1348 DO jmatrix = 1, mmatrix
1349 DO imatrix = 1, nmatrix
1350 NULLIFY (matrix_set(imatrix, jmatrix, kmatrix)%matrix)
1351 END DO
1352 END DO
1353 END DO
1354 END SUBROUTINE allocate_dbcsr_matrix_set_3d
1355
1356! **************************************************************************************************
1357!> \brief Allocate and initialize a real matrix 4-dimensional set.
1358!> \param[in,out] matrix_set Set containing the DBCSR matrix pointer type
1359!> \param[in] nmatrix Size of set
1360!> \param mmatrix ...
1361!> \param pmatrix ...
1362!> \param qmatrix ...
1363!> \par History
1364!> 2009-08-17 Adapted from sparse_matrix_type for DBCSR
1365! **************************************************************************************************
1366 SUBROUTINE allocate_dbcsr_matrix_set_4d(matrix_set, nmatrix, mmatrix, pmatrix, qmatrix)
1367 TYPE(dbcsr_p_type), DIMENSION(:, :, :, :), POINTER :: matrix_set
1368 INTEGER, INTENT(IN) :: nmatrix, mmatrix, pmatrix, qmatrix
1369
1370 INTEGER :: imatrix, jmatrix, kmatrix, lmatrix
1371
1372 IF (ASSOCIATED(matrix_set)) CALL dbcsr_deallocate_matrix_set(matrix_set)
1373 ALLOCATE (matrix_set(nmatrix, mmatrix, pmatrix, qmatrix))
1374 DO lmatrix = 1, qmatrix
1375 DO kmatrix = 1, pmatrix
1376 DO jmatrix = 1, mmatrix
1377 DO imatrix = 1, nmatrix
1378 NULLIFY (matrix_set(imatrix, jmatrix, kmatrix, lmatrix)%matrix)
1379 END DO
1380 END DO
1381 END DO
1382 END DO
1383 END SUBROUTINE allocate_dbcsr_matrix_set_4d
1384
1385! **************************************************************************************************
1386!> \brief Allocate and initialize a real matrix 5-dimensional set.
1387!> \param[in,out] matrix_set Set containing the DBCSR matrix pointer type
1388!> \param[in] nmatrix Size of set
1389!> \param mmatrix ...
1390!> \param pmatrix ...
1391!> \param qmatrix ...
1392!> \param smatrix ...
1393!> \par History
1394!> 2009-08-17 Adapted from sparse_matrix_type for DBCSR
1395! **************************************************************************************************
1396 SUBROUTINE allocate_dbcsr_matrix_set_5d(matrix_set, nmatrix, mmatrix, pmatrix, qmatrix, smatrix)
1397 TYPE(dbcsr_p_type), DIMENSION(:, :, :, :, :), &
1398 POINTER :: matrix_set
1399 INTEGER, INTENT(IN) :: nmatrix, mmatrix, pmatrix, qmatrix, &
1400 smatrix
1401
1402 INTEGER :: hmatrix, imatrix, jmatrix, kmatrix, &
1403 lmatrix
1404
1405 IF (ASSOCIATED(matrix_set)) CALL dbcsr_deallocate_matrix_set(matrix_set)
1406 ALLOCATE (matrix_set(nmatrix, mmatrix, pmatrix, qmatrix, smatrix))
1407 DO hmatrix = 1, smatrix
1408 DO lmatrix = 1, qmatrix
1409 DO kmatrix = 1, pmatrix
1410 DO jmatrix = 1, mmatrix
1411 DO imatrix = 1, nmatrix
1412 NULLIFY (matrix_set(imatrix, jmatrix, kmatrix, lmatrix, hmatrix)%matrix)
1413 END DO
1414 END DO
1415 END DO
1416 END DO
1417 END DO
1418 END SUBROUTINE allocate_dbcsr_matrix_set_5d
1419
1420 ! **************************************************************************************************
1421!> \brief Deallocate a real matrix set and release all of the member matrices.
1422!> \param[in,out] matrix_set Set containing the DBCSR matrix pointer type
1423!> \par History
1424!> 2009-08-17 Adapted from sparse_matrix_type for DBCSR
1425! **************************************************************************************************
1426 SUBROUTINE deallocate_dbcsr_matrix_set_1d(matrix_set)
1427
1428 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_set
1429
1430 INTEGER :: imatrix
1431
1432 IF (ASSOCIATED(matrix_set)) THEN
1433 DO imatrix = 1, SIZE(matrix_set)
1434 CALL dbcsr_deallocate_matrix(matrix_set(imatrix)%matrix)
1435 END DO
1436 DEALLOCATE (matrix_set)
1437 END IF
1438
1439 END SUBROUTINE deallocate_dbcsr_matrix_set_1d
1440
1441! **************************************************************************************************
1442!> \brief Deallocate a real matrix set and release all of the member matrices.
1443!> \param[in,out] matrix_set Set containing the DBCSR matrix pointer type
1444!> \par History
1445!> 2009-08-17 Adapted from sparse_matrix_type for DBCSR
1446! **************************************************************************************************
1447 SUBROUTINE deallocate_dbcsr_matrix_set_2d(matrix_set)
1448
1449 TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: matrix_set
1450
1451 INTEGER :: imatrix, jmatrix
1452
1453 IF (ASSOCIATED(matrix_set)) THEN
1454 DO jmatrix = 1, SIZE(matrix_set, 2)
1455 DO imatrix = 1, SIZE(matrix_set, 1)
1456 CALL dbcsr_deallocate_matrix(matrix_set(imatrix, jmatrix)%matrix)
1457 END DO
1458 END DO
1459 DEALLOCATE (matrix_set)
1460 END IF
1461 END SUBROUTINE deallocate_dbcsr_matrix_set_2d
1462
1463! **************************************************************************************************
1464!> \brief Deallocate a real matrix set and release all of the member matrices.
1465!> \param[in,out] matrix_set Set containing the DBCSR matrix pointer type
1466!> \par History
1467!> 2009-08-17 Adapted from sparse_matrix_type for DBCSR
1468! **************************************************************************************************
1469 SUBROUTINE deallocate_dbcsr_matrix_set_3d(matrix_set)
1470
1471 TYPE(dbcsr_p_type), DIMENSION(:, :, :), POINTER :: matrix_set
1472
1473 INTEGER :: imatrix, jmatrix, kmatrix
1474
1475 IF (ASSOCIATED(matrix_set)) THEN
1476 DO kmatrix = 1, SIZE(matrix_set, 3)
1477 DO jmatrix = 1, SIZE(matrix_set, 2)
1478 DO imatrix = 1, SIZE(matrix_set, 1)
1479 CALL dbcsr_deallocate_matrix(matrix_set(imatrix, jmatrix, kmatrix)%matrix)
1480 END DO
1481 END DO
1482 END DO
1483 DEALLOCATE (matrix_set)
1484 END IF
1485 END SUBROUTINE deallocate_dbcsr_matrix_set_3d
1486
1487! **************************************************************************************************
1488!> \brief Deallocate a real matrix set and release all of the member matrices.
1489!> \param[in,out] matrix_set Set containing the DBCSR matrix pointer type
1490!> \par History
1491!> 2009-08-17 Adapted from sparse_matrix_type for DBCSR
1492! **************************************************************************************************
1493 SUBROUTINE deallocate_dbcsr_matrix_set_4d(matrix_set)
1494
1495 TYPE(dbcsr_p_type), DIMENSION(:, :, :, :), POINTER :: matrix_set
1496
1497 INTEGER :: imatrix, jmatrix, kmatrix, lmatrix
1498
1499 IF (ASSOCIATED(matrix_set)) THEN
1500 DO lmatrix = 1, SIZE(matrix_set, 4)
1501 DO kmatrix = 1, SIZE(matrix_set, 3)
1502 DO jmatrix = 1, SIZE(matrix_set, 2)
1503 DO imatrix = 1, SIZE(matrix_set, 1)
1504 CALL dbcsr_deallocate_matrix(matrix_set(imatrix, jmatrix, kmatrix, lmatrix)%matrix)
1505 END DO
1506 END DO
1507 END DO
1508 END DO
1509 DEALLOCATE (matrix_set)
1510 END IF
1511 END SUBROUTINE deallocate_dbcsr_matrix_set_4d
1512
1513! **************************************************************************************************
1514!> \brief Deallocate a real matrix set and release all of the member matrices.
1515!> \param[in,out] matrix_set Set containing the DBCSR matrix pointer type
1516!> \par History
1517!> 2009-08-17 Adapted from sparse_matrix_type for DBCSR
1518! **************************************************************************************************
1519 SUBROUTINE deallocate_dbcsr_matrix_set_5d(matrix_set)
1520
1521 TYPE(dbcsr_p_type), DIMENSION(:, :, :, :, :), &
1522 POINTER :: matrix_set
1523
1524 INTEGER :: hmatrix, imatrix, jmatrix, kmatrix, &
1525 lmatrix
1526
1527 IF (ASSOCIATED(matrix_set)) THEN
1528 DO hmatrix = 1, SIZE(matrix_set, 5)
1529 DO lmatrix = 1, SIZE(matrix_set, 4)
1530 DO kmatrix = 1, SIZE(matrix_set, 3)
1531 DO jmatrix = 1, SIZE(matrix_set, 2)
1532 DO imatrix = 1, SIZE(matrix_set, 1)
1533 CALL dbcsr_deallocate_matrix(matrix_set(imatrix, jmatrix, kmatrix, lmatrix, hmatrix)%matrix)
1534 END DO
1535 END DO
1536 END DO
1537 END DO
1538 END DO
1539 DEALLOCATE (matrix_set)
1540 END IF
1541 END SUBROUTINE deallocate_dbcsr_matrix_set_5d
1542
1543END MODULE cp_dbcsr_operations
dgemm
static void dgemm(const char transa, const char transb, const int m, const int n, const int k, const double alpha, const double *a, const int lda, const double *b, const int ldb, const double beta, double *c, const int ldc)
Convenient wrapper to hide Fortran nature of dgemm_, swapping a and b.
Definition grid_cpu_task_list.c:206
cp_dbcsr_api::dbcsr_create
Definition cp_dbcsr_api.F:192
cp_dbcsr_operations::dbcsr_allocate_matrix_set
Definition cp_dbcsr_operations.F:77
cp_dbcsr_operations::dbcsr_deallocate_matrix_set
Definition cp_dbcsr_operations.F:85
cp_fm_types::cp_fm_release
Definition cp_fm_types.F:87
cp_fm_types::cp_fm_to_fm
Definition cp_fm_types.F:82
cp_blacs_env
methods related to the blacs parallel environment
Definition cp_blacs_env.F:15
cp_dbcsr_api
Definition cp_dbcsr_api.F:8
cp_dbcsr_api::dbcsr_verify_matrix
subroutine, public dbcsr_verify_matrix(matrix, verbosity, local)
...
Definition cp_dbcsr_api.F:1264
cp_dbcsr_api::dbcsr_distribution_release
subroutine, public dbcsr_distribution_release(dist)
...
Definition cp_dbcsr_api.F:629
cp_dbcsr_api::dbcsr_scale
subroutine, public dbcsr_scale(matrix, alpha_scalar)
...
Definition cp_dbcsr_api.F:1165
cp_dbcsr_api::dbcsr_distribution_new
subroutine, public dbcsr_distribution_new(dist, template, group, pgrid, row_dist, col_dist, reuse_arrays)
...
Definition cp_dbcsr_api.F:608
cp_dbcsr_api::dbcsr_deallocate_matrix
subroutine, public dbcsr_deallocate_matrix(matrix)
...
Definition cp_dbcsr_api.F:233
cp_dbcsr_api::dbcsr_iterator_next_block
subroutine, public dbcsr_iterator_next_block(iterator, row, column, block, block_number_argument_has_been_removed, row_size, col_size, row_offset, col_offset)
...
Definition cp_dbcsr_api.F:956
cp_dbcsr_api::dbcsr_get_matrix_type
character function, public dbcsr_get_matrix_type(matrix)
...
Definition cp_dbcsr_api.F:847
cp_dbcsr_api::dbcsr_iterator_blocks_left
logical function, public dbcsr_iterator_blocks_left(iterator)
...
Definition cp_dbcsr_api.F:930
cp_dbcsr_api::dbcsr_iterator_stop
subroutine, public dbcsr_iterator_stop(iterator)
...
Definition cp_dbcsr_api.F:1027
cp_dbcsr_api::dbcsr_valid_index
logical function, public dbcsr_valid_index(matrix)
...
Definition cp_dbcsr_api.F:1247
cp_dbcsr_api::dbcsr_desymmetrize
subroutine, public dbcsr_desymmetrize(matrix_a, matrix_b)
...
Definition cp_dbcsr_api.F:517
cp_dbcsr_api::dbcsr_copy
subroutine, public dbcsr_copy(matrix_b, matrix_a, name, keep_sparsity, keep_imaginary)
...
Definition cp_dbcsr_api.F:368
cp_dbcsr_api::dbcsr_multiply
subroutine, public dbcsr_multiply(transa, transb, alpha, matrix_a, matrix_b, beta, matrix_c, first_row, last_row, first_column, last_column, first_k, last_k, retain_sparsity, filter_eps, flop)
...
Definition cp_dbcsr_api.F:1073
cp_dbcsr_api::dbcsr_get_info
subroutine, public dbcsr_get_info(matrix, nblkrows_total, nblkcols_total, nfullrows_total, nfullcols_total, nblkrows_local, nblkcols_local, nfullrows_local, nfullcols_local, my_prow, my_pcol, local_rows, local_cols, proc_row_dist, proc_col_dist, row_blk_size, col_blk_size, row_blk_offset, col_blk_offset, distribution, name, matrix_type, group)
...
Definition cp_dbcsr_api.F:794
cp_dbcsr_api::dbcsr_iterator_start
subroutine, public dbcsr_iterator_start(iterator, matrix, shared, dynamic, dynamic_byrows)
...
Definition cp_dbcsr_api.F:989
cp_dbcsr_api::dbcsr_release
subroutine, public dbcsr_release(matrix)
...
Definition cp_dbcsr_api.F:1119
cp_dbcsr_api::dbcsr_complete_redistribute
subroutine, public dbcsr_complete_redistribute(matrix, redist)
...
Definition cp_dbcsr_api.F:317
cp_dbcsr_api::dbcsr_iterator_readonly_start
subroutine, public dbcsr_iterator_readonly_start(iterator, matrix, shared, dynamic, dynamic_byrows)
Like dbcsr_iterator_start() but with matrix being INTENT(IN). When invoking this routine,...
Definition cp_dbcsr_api.F:1010
cp_dbcsr_api::dbcsr_add
subroutine, public dbcsr_add(matrix_a, matrix_b, alpha_scalar, beta_scalar)
...
Definition cp_dbcsr_api.F:251
cp_dbcsr_api::dbcsr_distribution_get
subroutine, public dbcsr_distribution_get(dist, row_dist, col_dist, nrows, ncols, has_threads, group, mynode, numnodes, nprows, npcols, myprow, mypcol, pgrid, subgroups_defined, prow_group, pcol_group)
...
Definition cp_dbcsr_api.F:564
cp_dbcsr_contrib
Definition cp_dbcsr_contrib.F:8
cp_dbcsr_contrib::dbcsr_reserve_all_blocks
subroutine, public dbcsr_reserve_all_blocks(matrix)
Reserves all blocks.
Definition cp_dbcsr_contrib.F:295
cp_dbcsr_contrib::dbcsr_frobenius_norm
real(dp) function, public dbcsr_frobenius_norm(matrix)
Compute the frobenius norm of a dbcsr matrix.
Definition cp_dbcsr_contrib.F:127
cp_dbcsr_operations
DBCSR operations in CP2K.
Definition cp_dbcsr_operations.F:18
cp_dbcsr_operations::cp_dbcsr_sm_fm_multiply
subroutine, public cp_dbcsr_sm_fm_multiply(matrix, fm_in, fm_out, ncol, alpha, beta)
multiply a dbcsr with a fm matrix
Definition cp_dbcsr_operations.F:552
cp_dbcsr_operations::max_elements_per_block
integer, save, public max_elements_per_block
Definition cp_dbcsr_operations.F:57
cp_dbcsr_operations::dbcsr_create_dist_r_unrot
subroutine, public dbcsr_create_dist_r_unrot(dist_right, dist_left, ncolumns, right_col_blk_sizes)
Creates a new distribution for the right matrix in a matrix multiplication with unrotated grid.
Definition cp_dbcsr_operations.F:1098
cp_dbcsr_operations::dbcsr_multiply_local
subroutine, public dbcsr_multiply_local(matrix_a, vec_b, vec_c, ncol, alpha)
multiply a dbcsr with a replicated array c = alpha_scalar * A (dbscr) * b + c
Definition cp_dbcsr_operations.F:459
cp_dbcsr_operations::cp_dbcsr_dist2d_to_dist
subroutine, public cp_dbcsr_dist2d_to_dist(dist2d, dist)
Creates a DBCSR distribution from a distribution_2d.
Definition cp_dbcsr_operations.F:419
cp_dbcsr_operations::copy_dbcsr_to_fm
subroutine, public copy_dbcsr_to_fm(matrix, fm)
Copy a DBCSR matrix to a BLACS matrix.
Definition cp_dbcsr_operations.F:214
cp_dbcsr_operations::copy_fm_to_dbcsr_bc
subroutine, public copy_fm_to_dbcsr_bc(fm, bc_mat)
Copy a BLACS matrix to a dbcsr matrix with a special block-cyclic distribution, which requires no com...
Definition cp_dbcsr_operations.F:149
cp_dbcsr_operations::copy_dbcsr_to_fm_bc
subroutine, public copy_dbcsr_to_fm_bc(bc_mat, fm)
Copy a DBCSR_BLACS matrix to a BLACS matrix.
Definition cp_dbcsr_operations.F:275
cp_dbcsr_operations::cp_fm_to_dbcsr_row_template
subroutine, public cp_fm_to_dbcsr_row_template(matrix, fm_in, template)
Utility function to copy a specially shaped fm to dbcsr_matrix The result matrix will be the matrix i...
Definition cp_dbcsr_operations.F:878
cp_dbcsr_operations::cp_dbcsr_m_by_n_from_row_template
subroutine, public cp_dbcsr_m_by_n_from_row_template(matrix, template, n, sym)
Utility function to create dbcsr matrix, m x n matrix (n arbitrary) with the same processor grid and ...
Definition cp_dbcsr_operations.F:955
cp_dbcsr_operations::cp_dbcsr_plus_fm_fm_t
subroutine, public cp_dbcsr_plus_fm_fm_t(sparse_matrix, matrix_v, matrix_g, ncol, alpha, keep_sparsity, symmetry_mode)
performs the multiplication sparse_matrix+dense_mat*dens_mat^T if matrix_g is not explicitly given,...
Definition cp_dbcsr_operations.F:695
cp_dbcsr_operations::cp_dbcsr_m_by_n_from_template
subroutine, public cp_dbcsr_m_by_n_from_template(matrix, template, m, n, sym)
Utility function to create an arbitrary shaped dbcsr matrix with the same processor grid as the templ...
Definition cp_dbcsr_operations.F:912
cp_dbcsr_operations::dbcsr_copy_columns_hack
subroutine, public dbcsr_copy_columns_hack(matrix_b, matrix_a, ncol, source_start, target_start, para_env, blacs_env)
hack for dbcsr_copy_columns
Definition cp_dbcsr_operations.F:375
cp_dbcsr_operations::copy_fm_to_dbcsr
subroutine, public copy_fm_to_dbcsr(fm, matrix, keep_sparsity)
Copy a BLACS matrix to a dbcsr matrix.
Definition cp_dbcsr_operations.F:111
cp_fm_basic_linalg
basic linear algebra operations for full matrices
Definition cp_fm_basic_linalg.F:14
cp_fm_basic_linalg::cp_fm_gemm
subroutine, public cp_fm_gemm(transa, transb, m, n, k, alpha, matrix_a, matrix_b, beta, matrix_c, a_first_col, a_first_row, b_first_col, b_first_row, c_first_col, c_first_row)
computes matrix_c = beta * matrix_c + alpha * ( matrix_a ** transa ) * ( matrix_b ** transb )
Definition cp_fm_basic_linalg.F:440
cp_fm_struct
represent the structure of a full matrix
Definition cp_fm_struct.F:14
cp_fm_struct::cp_fm_struct_create
subroutine, public cp_fm_struct_create(fmstruct, para_env, context, nrow_global, ncol_global, nrow_block, ncol_block, descriptor, first_p_pos, local_leading_dimension, template_fmstruct, square_blocks, force_block)
allocates and initializes a full matrix structure
Definition cp_fm_struct.F:132
cp_fm_struct::cp_fm_struct_release
subroutine, public cp_fm_struct_release(fmstruct)
releases a full matrix structure
Definition cp_fm_struct.F:351
cp_fm_types
represent a full matrix distributed on many processors
Definition cp_fm_types.F:15
cp_fm_types::cp_fm_get_info
subroutine, public cp_fm_get_info(matrix, name, nrow_global, ncol_global, nrow_block, ncol_block, nrow_local, ncol_local, row_indices, col_indices, local_data, context, nrow_locals, ncol_locals, matrix_struct, para_env)
returns all kind of information about the full matrix
Definition cp_fm_types.F:1009
cp_fm_types::cp_fm_create
subroutine, public cp_fm_create(matrix, matrix_struct, name, use_sp)
creates a new full matrix with the given structure
Definition cp_fm_types.F:164
distribution_2d_types
stores a mapping of 2D info (e.g. matrix) on a 2D processor distribution (i.e. blacs grid) where cpus...
Definition distribution_2d_types.F:15
distribution_2d_types::distribution_2d_get
subroutine, public distribution_2d_get(distribution_2d, row_distribution, col_distribution, n_row_distribution, n_col_distribution, n_local_rows, n_local_cols, local_rows, local_cols, flat_local_rows, flat_local_cols, n_flat_local_rows, n_flat_local_cols, blacs_env)
returns various attributes about the distribution_2d
Definition distribution_2d_types.F:481
kinds
Defines the basic variable types.
Definition kinds.F:23
kinds::dp
integer, parameter, public dp
Definition kinds.F:34
kinds::default_string_length
integer, parameter, public default_string_length
Definition kinds.F:57
mathlib
Collection of simple mathematical functions and subroutines.
Definition mathlib.F:15
mathlib::lcm
elemental integer function, public lcm(a, b)
computes the least common multiplier of two numbers
Definition mathlib.F:1315
mathlib::gcd
elemental integer function, public gcd(a, b)
computes the greatest common divisor of two number
Definition mathlib.F:1280
message_passing
Interface to the message passing library MPI.
Definition message_passing.F:23
cp_blacs_env::cp_blacs_env_type
represent a blacs multidimensional parallel environment (for the mpi corrispective see cp_paratypes/m...
Definition cp_blacs_env.F:53
cp_dbcsr_api::dbcsr_distribution_type
Definition cp_dbcsr_api.F:180
cp_dbcsr_api::dbcsr_iterator_type
Definition cp_dbcsr_api.F:186
cp_dbcsr_api::dbcsr_p_type
Definition cp_dbcsr_api.F:170
cp_dbcsr_api::dbcsr_type
Definition cp_dbcsr_api.F:174
cp_fm_struct::cp_fm_struct_type
keeps the information about the structure of a full matrix
Definition cp_fm_struct.F:82
cp_fm_types::cp_fm_type
represent a full matrix
Definition cp_fm_types.F:113
distribution_2d_types::distribution_2d_type
distributes pairs on a 2d grid of processors
Definition distribution_2d_types.F:62
message_passing::mp_para_env_type
stores all the informations relevant to an mpi environment
Definition message_passing.F:763