22 USE dbcsr_api,
ONLY: dbcsr_distribution_get, &
23 dbcsr_convert_sizes_to_offsets, dbcsr_add, &
24 dbcsr_complete_redistribute, dbcsr_copy, dbcsr_create, &
25 dbcsr_deallocate_matrix, &
26 dbcsr_desymmetrize, dbcsr_distribution_new, &
27 dbcsr_get_data_type, dbcsr_get_info, dbcsr_get_matrix_type, &
28 dbcsr_iterator_type, dbcsr_iterator_blocks_left, dbcsr_iterator_next_block, &
29 dbcsr_iterator_start, dbcsr_iterator_stop, &
30 dbcsr_multiply, dbcsr_norm, dbcsr_p_type, dbcsr_release, &
31 dbcsr_reserve_all_blocks, dbcsr_scale, dbcsr_type, &
32 dbcsr_valid_index, dbcsr_verify_matrix, &
33 dbcsr_distribution_type, dbcsr_distribution_release, &
34 dbcsr_norm_frobenius, &
35 dbcsr_type_antisymmetric, dbcsr_type_complex_8, dbcsr_type_no_symmetry, dbcsr_type_real_8, &
53 #include "base/base_uses.f90"
57 CHARACTER(len=*),
PARAMETER,
PRIVATE :: moduleN =
'cp_dbcsr_operations'
82 MODULE PROCEDURE allocate_dbcsr_matrix_set_1d
83 MODULE PROCEDURE allocate_dbcsr_matrix_set_2d
84 MODULE PROCEDURE allocate_dbcsr_matrix_set_3d
85 MODULE PROCEDURE allocate_dbcsr_matrix_set_4d
86 MODULE PROCEDURE allocate_dbcsr_matrix_set_5d
90 MODULE PROCEDURE deallocate_dbcsr_matrix_set_1d
91 MODULE PROCEDURE deallocate_dbcsr_matrix_set_2d
92 MODULE PROCEDURE deallocate_dbcsr_matrix_set_3d
93 MODULE PROCEDURE deallocate_dbcsr_matrix_set_4d
94 MODULE PROCEDURE deallocate_dbcsr_matrix_set_5d
118 TYPE(cp_fm_type),
INTENT(IN) :: fm
119 TYPE(dbcsr_type),
INTENT(INOUT) :: matrix
120 LOGICAL,
INTENT(IN),
OPTIONAL :: keep_sparsity
122 CHARACTER(LEN=*),
PARAMETER :: routinen =
'copy_fm_to_dbcsr'
124 TYPE(dbcsr_type) :: bc_mat
127 CALL timeset(routinen, handle)
130 CALL dbcsr_complete_redistribute(bc_mat, matrix, keep_sparsity=keep_sparsity)
131 CALL dbcsr_release(bc_mat)
133 CALL timestop(handle)
143 TYPE(cp_fm_type),
INTENT(IN) :: fm
144 TYPE(dbcsr_type),
INTENT(INOUT) :: bc_mat
146 CHARACTER(LEN=*),
PARAMETER :: routinen =
'copy_fm_to_dbcsr_bc'
148 INTEGER :: col, handle, ncol_block, ncol_global, nrow_block, nrow_global, row
149 INTEGER,
ALLOCATABLE,
DIMENSION(:) :: first_col, first_row, last_col, last_row
150 INTEGER,
DIMENSION(:),
POINTER :: col_blk_size, row_blk_size
151 REAL (kind=
dp),
DIMENSION(:, :),
POINTER :: fm_block, dbcsr_block
152 TYPE(dbcsr_distribution_type) :: bc_dist
153 TYPE(dbcsr_iterator_type) :: iter
154 INTEGER,
DIMENSION(:, :),
POINTER :: pgrid
156 CALL timeset(routinen, handle)
158 IF (fm%use_sp) cpabort(
"copy_fm_to_dbcsr_bc: single precision not supported")
161 pgrid => fm%matrix_struct%context%blacs2mpi
164 nrow_block = fm%matrix_struct%nrow_block
165 ncol_block = fm%matrix_struct%ncol_block
166 nrow_global = fm%matrix_struct%nrow_global
167 ncol_global = fm%matrix_struct%ncol_global
168 NULLIFY (col_blk_size, row_blk_size)
169 CALL dbcsr_create_dist_block_cyclic(bc_dist, &
170 nrows=nrow_global, ncolumns=ncol_global, &
171 nrow_block=nrow_block, ncol_block=ncol_block, &
172 group_handle=fm%matrix_struct%para_env%get_handle(), pgrid=pgrid, &
173 row_blk_sizes=row_blk_size, col_blk_sizes=col_blk_size)
176 CALL dbcsr_create(bc_mat,
"Block-cyclic ", bc_dist, &
177 dbcsr_type_no_symmetry, row_blk_size, col_blk_size, nze=0, &
178 reuse_arrays=.true., data_type=dbcsr_type_real_8)
179 CALL dbcsr_distribution_release(bc_dist)
182 CALL dbcsr_reserve_all_blocks(bc_mat)
184 CALL calculate_fm_block_ranges(bc_mat, first_row, last_row, first_col, last_col)
188 fm_block => fm%local_data
191 CALL dbcsr_iterator_start(iter, bc_mat)
192 DO WHILE (dbcsr_iterator_blocks_left(iter))
193 CALL dbcsr_iterator_next_block(iter, row, col, dbcsr_block)
194 dbcsr_block(:, :) = fm_block(first_row(row):last_row(row), first_col(col):last_col(col))
196 CALL dbcsr_iterator_stop(iter)
199 CALL timestop(handle)
208 TYPE(dbcsr_type),
INTENT(IN) :: matrix
209 TYPE(cp_fm_type),
INTENT(IN) :: fm
211 CHARACTER(LEN=*),
PARAMETER :: routinen =
'copy_dbcsr_to_fm'
213 INTEGER,
DIMENSION(:),
POINTER :: col_blk_size, row_blk_size
214 INTEGER :: handle, ncol_block, nfullcols_total, group_handle, &
215 nfullrows_total, nrow_block
216 TYPE(dbcsr_type) :: bc_mat, matrix_nosym
217 TYPE(dbcsr_distribution_type) :: dist, bc_dist
218 CHARACTER(len=default_string_length) :: name
219 INTEGER,
DIMENSION(:, :),
POINTER :: pgrid
221 CALL timeset(routinen, handle)
224 CALL dbcsr_get_info(matrix, &
227 nfullrows_total=nfullrows_total, &
228 nfullcols_total=nfullcols_total)
230 cpassert(fm%matrix_struct%nrow_global == nfullrows_total)
231 cpassert(fm%matrix_struct%ncol_global == nfullcols_total)
234 nrow_block = fm%matrix_struct%nrow_block
235 ncol_block = fm%matrix_struct%ncol_block
238 NULLIFY (col_blk_size, row_blk_size)
239 CALL dbcsr_distribution_get(dist, group=group_handle, pgrid=pgrid)
240 CALL dbcsr_create_dist_block_cyclic(bc_dist, &
241 nrows=nfullrows_total, ncolumns=nfullcols_total, &
242 nrow_block=nrow_block, ncol_block=ncol_block, &
243 group_handle=group_handle, pgrid=pgrid, &
244 row_blk_sizes=row_blk_size, col_blk_sizes=col_blk_size)
246 CALL dbcsr_create(bc_mat,
"Block-cyclic"//name, bc_dist, &
247 dbcsr_type_no_symmetry, row_blk_size, col_blk_size, &
248 nze=0, data_type=dbcsr_get_data_type(matrix), &
250 CALL dbcsr_distribution_release(bc_dist)
252 CALL dbcsr_create(matrix_nosym, template=matrix, matrix_type=
"N")
253 CALL dbcsr_desymmetrize(matrix, matrix_nosym)
254 CALL dbcsr_complete_redistribute(matrix_nosym, bc_mat)
255 CALL dbcsr_release(matrix_nosym)
259 CALL dbcsr_release(bc_mat)
261 CALL timestop(handle)
270 TYPE(dbcsr_type),
INTENT(IN) :: bc_mat
271 TYPE(cp_fm_type),
INTENT(IN) :: fm
273 CHARACTER(LEN=*),
PARAMETER :: routinen =
'copy_dbcsr_to_fm_bc'
275 INTEGER :: col, handle, row
276 INTEGER,
ALLOCATABLE,
DIMENSION(:) :: first_col, first_row, last_col, last_row
277 REAL (kind=
dp),
DIMENSION(:, :),
POINTER :: dbcsr_block, fm_block
278 TYPE(dbcsr_iterator_type) :: iter
280 CALL timeset(routinen, handle)
282 IF (fm%use_sp) cpabort(
"copy_dbcsr_to_fm_bc: single precision not supported")
284 CALL calculate_fm_block_ranges(bc_mat, first_row, last_row, first_col, last_col)
287 fm_block => fm%local_data
288 fm_block = real(0.0, kind=
dp)
291 CALL dbcsr_iterator_start(iter, bc_mat)
292 DO WHILE (dbcsr_iterator_blocks_left(iter))
293 CALL dbcsr_iterator_next_block(iter, row, col, dbcsr_block)
294 fm_block(first_row(row):last_row(row), first_col(col):last_col(col)) = dbcsr_block(:, :)
296 CALL dbcsr_iterator_stop(iter)
299 CALL timestop(handle)
319 TYPE(cp_cfm_type),
INTENT(IN) :: fm
320 TYPE(dbcsr_type),
INTENT(INOUT) :: matrix
321 LOGICAL,
INTENT(IN),
OPTIONAL :: keep_sparsity
323 CHARACTER(LEN=*),
PARAMETER :: routinen =
'copy_cfm_to_dbcsr'
325 TYPE(dbcsr_type) :: bc_mat
328 CALL timeset(routinen, handle)
330 CALL copy_cfm_to_dbcsr_bc(fm, bc_mat)
331 CALL dbcsr_complete_redistribute(bc_mat, matrix, keep_sparsity=keep_sparsity)
332 CALL dbcsr_release(bc_mat)
334 CALL timestop(handle)
343 SUBROUTINE copy_cfm_to_dbcsr_bc(fm, bc_mat)
344 TYPE(cp_cfm_type),
INTENT(IN) :: fm
345 TYPE(dbcsr_type),
INTENT(INOUT) :: bc_mat
347 CHARACTER(LEN=*),
PARAMETER :: routinen =
'copy_cfm_to_dbcsr_bc'
349 INTEGER :: col, handle, ncol_block, ncol_global, nrow_block, nrow_global, row
350 INTEGER,
ALLOCATABLE,
DIMENSION(:) :: first_col, first_row, last_col, last_row
351 INTEGER,
DIMENSION(:),
POINTER :: col_blk_size, row_blk_size
352 COMPLEX (KIND=dp),
DIMENSION(:, :),
POINTER :: fm_block, dbcsr_block
353 TYPE(dbcsr_distribution_type) :: bc_dist
354 TYPE(dbcsr_iterator_type) :: iter
355 INTEGER,
DIMENSION(:, :),
POINTER :: pgrid
357 CALL timeset(routinen, handle)
361 pgrid => fm%matrix_struct%context%blacs2mpi
364 nrow_block = fm%matrix_struct%nrow_block
365 ncol_block = fm%matrix_struct%ncol_block
366 nrow_global = fm%matrix_struct%nrow_global
367 ncol_global = fm%matrix_struct%ncol_global
368 NULLIFY (col_blk_size, row_blk_size)
369 CALL dbcsr_create_dist_block_cyclic(bc_dist, &
370 nrows=nrow_global, ncolumns=ncol_global, &
371 nrow_block=nrow_block, ncol_block=ncol_block, &
372 group_handle=fm%matrix_struct%para_env%get_handle(), pgrid=pgrid, &
373 row_blk_sizes=row_blk_size, col_blk_sizes=col_blk_size)
376 CALL dbcsr_create(bc_mat,
"Block-cyclic ", bc_dist, &
377 dbcsr_type_no_symmetry, row_blk_size, col_blk_size, nze=0, &
378 reuse_arrays=.true., data_type=dbcsr_type_complex_8)
379 CALL dbcsr_distribution_release(bc_dist)
382 CALL dbcsr_reserve_all_blocks(bc_mat)
384 CALL calculate_fm_block_ranges(bc_mat, first_row, last_row, first_col, last_col)
388 fm_block => fm%local_data
391 CALL dbcsr_iterator_start(iter, bc_mat)
392 DO WHILE (dbcsr_iterator_blocks_left(iter))
393 CALL dbcsr_iterator_next_block(iter, row, col, dbcsr_block)
394 dbcsr_block(:, :) = fm_block(first_row(row):last_row(row), first_col(col):last_col(col))
396 CALL dbcsr_iterator_stop(iter)
399 CALL timestop(handle)
400 END SUBROUTINE copy_cfm_to_dbcsr_bc
408 TYPE(dbcsr_type),
INTENT(IN) :: matrix
409 TYPE(cp_cfm_type),
INTENT(IN) :: fm
411 CHARACTER(LEN=*),
PARAMETER :: routinen =
'copy_dbcsr_to_cfm'
413 INTEGER,
DIMENSION(:),
POINTER :: col_blk_size, row_blk_size
414 INTEGER :: handle, ncol_block, nfullcols_total, group_handle, &
415 nfullrows_total, nrow_block
416 TYPE(dbcsr_type) :: bc_mat, matrix_nosym
417 TYPE(dbcsr_distribution_type) :: dist, bc_dist
418 CHARACTER(len=default_string_length) :: name
419 INTEGER,
DIMENSION(:, :),
POINTER :: pgrid
421 CALL timeset(routinen, handle)
424 CALL dbcsr_get_info(matrix, &
427 nfullrows_total=nfullrows_total, &
428 nfullcols_total=nfullcols_total)
430 cpassert(fm%matrix_struct%nrow_global == nfullrows_total)
431 cpassert(fm%matrix_struct%ncol_global == nfullcols_total)
434 nrow_block = fm%matrix_struct%nrow_block
435 ncol_block = fm%matrix_struct%ncol_block
438 NULLIFY (col_blk_size, row_blk_size)
439 CALL dbcsr_distribution_get(dist, group=group_handle, pgrid=pgrid)
440 CALL dbcsr_create_dist_block_cyclic(bc_dist, &
441 nrows=nfullrows_total, ncolumns=nfullcols_total, &
442 nrow_block=nrow_block, ncol_block=ncol_block, &
443 group_handle=group_handle, pgrid=pgrid, &
444 row_blk_sizes=row_blk_size, col_blk_sizes=col_blk_size)
446 CALL dbcsr_create(bc_mat,
"Block-cyclic"//name, bc_dist, &
447 dbcsr_type_no_symmetry, row_blk_size, col_blk_size, &
448 nze=0, data_type=dbcsr_get_data_type(matrix), &
450 CALL dbcsr_distribution_release(bc_dist)
452 CALL dbcsr_create(matrix_nosym, template=matrix, matrix_type=
"N")
453 CALL dbcsr_desymmetrize(matrix, matrix_nosym)
454 CALL dbcsr_complete_redistribute(matrix_nosym, bc_mat)
455 CALL dbcsr_release(matrix_nosym)
457 CALL copy_dbcsr_to_cfm_bc(bc_mat, fm)
459 CALL dbcsr_release(bc_mat)
461 CALL timestop(handle)
469 SUBROUTINE copy_dbcsr_to_cfm_bc(bc_mat, fm)
470 TYPE(dbcsr_type),
INTENT(IN) :: bc_mat
471 TYPE(cp_cfm_type),
INTENT(IN) :: fm
473 CHARACTER(LEN=*),
PARAMETER :: routinen =
'copy_dbcsr_to_cfm_bc'
475 INTEGER :: col, handle, row
476 INTEGER,
ALLOCATABLE,
DIMENSION(:) :: first_col, first_row, last_col, last_row
477 COMPLEX (KIND=dp),
DIMENSION(:, :),
POINTER :: dbcsr_block, fm_block
478 TYPE(dbcsr_iterator_type) :: iter
480 CALL timeset(routinen, handle)
483 CALL calculate_fm_block_ranges(bc_mat, first_row, last_row, first_col, last_col)
486 fm_block => fm%local_data
487 fm_block = cmplx(0.0, kind=
dp)
490 CALL dbcsr_iterator_start(iter, bc_mat)
491 DO WHILE (dbcsr_iterator_blocks_left(iter))
492 CALL dbcsr_iterator_next_block(iter, row, col, dbcsr_block)
493 fm_block(first_row(row):last_row(row), first_col(col):last_col(col)) = dbcsr_block(:, :)
495 CALL dbcsr_iterator_stop(iter)
498 CALL timestop(handle)
499 END SUBROUTINE copy_dbcsr_to_cfm_bc
511 SUBROUTINE calculate_fm_block_ranges(bc_mat, first_row, last_row, first_col, last_col)
512 TYPE(dbcsr_type),
INTENT(IN) :: bc_mat
513 INTEGER :: col, nblkcols_local, nblkcols_total, nblkrows_local, nblkrows_total, row
514 INTEGER,
ALLOCATABLE,
DIMENSION(:) :: first_col, first_row, last_col, &
515 last_row, local_col_sizes, &
517 INTEGER,
DIMENSION(:),
POINTER :: col_blk_size, local_cols, local_rows, &
520 CALL dbcsr_get_info(bc_mat, &
521 nblkrows_total=nblkrows_total, &
522 nblkcols_total=nblkcols_total, &
523 nblkrows_local=nblkrows_local, &
524 nblkcols_local=nblkcols_local, &
525 local_rows=local_rows, &
526 local_cols=local_cols, &
527 row_blk_size=row_blk_size, &
528 col_blk_size=col_blk_size)
531 ALLOCATE (local_row_sizes(nblkrows_total))
532 local_row_sizes(:) = 0
533 IF (nblkrows_local .GE. 1)
THEN
534 DO row = 1, nblkrows_local
535 local_row_sizes(local_rows(row)) = row_blk_size(local_rows(row))
538 ALLOCATE (first_row(nblkrows_total), last_row(nblkrows_total))
539 CALL dbcsr_convert_sizes_to_offsets(local_row_sizes, first_row, last_row)
542 ALLOCATE (local_col_sizes(nblkcols_total))
543 local_col_sizes(:) = 0
544 IF (nblkcols_local .GE. 1)
THEN
545 DO col = 1, nblkcols_local
546 local_col_sizes(local_cols(col)) = col_blk_size(local_cols(col))
549 ALLOCATE (first_col(nblkcols_total), last_col(nblkcols_total))
550 CALL dbcsr_convert_sizes_to_offsets(local_col_sizes, first_col, last_col)
552 END SUBROUTINE calculate_fm_block_ranges
566 ncol, source_start, target_start, para_env, blacs_env)
568 TYPE(dbcsr_type),
INTENT(INOUT) :: matrix_b
569 TYPE(dbcsr_type),
INTENT(IN) :: matrix_a
570 INTEGER,
INTENT(IN) :: ncol, source_start, target_start
571 TYPE(mp_para_env_type),
POINTER :: para_env
572 TYPE(cp_blacs_env_type),
POINTER :: blacs_env
574 INTEGER :: nfullcols_total, nfullrows_total
575 TYPE(cp_fm_struct_type),
POINTER :: fm_struct
576 TYPE(cp_fm_type) :: fm_matrix_a, fm_matrix_b
579 CALL dbcsr_get_info(matrix_a, nfullrows_total=nfullrows_total, nfullcols_total=nfullcols_total)
581 ncol_global=nfullcols_total, para_env=para_env)
582 CALL cp_fm_create(fm_matrix_a, fm_struct, name=
"fm_matrix_a")
585 CALL dbcsr_get_info(matrix_b, nfullrows_total=nfullrows_total, nfullcols_total=nfullcols_total)
587 ncol_global=nfullcols_total, para_env=para_env)
588 CALL cp_fm_create(fm_matrix_b, fm_struct, name=
"fm_matrix_b")
594 CALL cp_fm_to_fm(fm_matrix_a, fm_matrix_b, ncol, source_start, target_start)
598 CALL cp_fm_release(fm_matrix_a)
599 CALL cp_fm_release(fm_matrix_b)
611 TYPE(distribution_2d_type),
INTENT(IN),
TARGET :: dist2d
612 TYPE(dbcsr_distribution_type),
INTENT(OUT) :: dist
614 INTEGER,
DIMENSION(:, :),
POINTER :: pgrid, col_dist_data, row_dist_data
615 TYPE(cp_blacs_env_type),
POINTER :: blacs_env
616 TYPE(mp_para_env_type),
POINTER :: para_env
617 TYPE(distribution_2d_type),
POINTER :: dist2d_p
618 INTEGER,
DIMENSION(:),
POINTER :: row_dist, col_dist
622 row_distribution=row_dist_data, &
623 col_distribution=col_dist_data, &
625 CALL blacs_env%get(para_env=para_env, blacs2mpi=pgrid)
628 row_dist => row_dist_data(:, 1)
629 col_dist => col_dist_data(:, 1)
633 CALL dbcsr_distribution_new(dist, &
634 group=para_env%get_handle(), pgrid=pgrid, &
651 TYPE(dbcsr_type),
INTENT(IN) :: matrix_a
652 REAL(
dp),
DIMENSION(:, :),
INTENT(IN) :: vec_b
653 REAL(
dp),
DIMENSION(:, :),
INTENT(INOUT) :: vec_c
654 INTEGER,
INTENT(in),
OPTIONAL :: ncol
655 REAL(
dp),
INTENT(IN),
OPTIONAL :: alpha
657 CHARACTER(LEN=*),
PARAMETER :: routinen =
'dbcsr_multiply_local'
659 INTEGER :: blk, col, coloff, my_ncol, row, rowoff, &
662 REAL(
dp) :: my_alpha, my_alpha2
663 REAL(
dp),
DIMENSION(:, :),
POINTER :: data_d
664 TYPE(dbcsr_iterator_type) :: iter
666 CALL timeset(routinen, timing_handle)
669 IF (
PRESENT(alpha)) my_alpha = alpha
671 my_ncol =
SIZE(vec_b, 2)
672 IF (
PRESENT(ncol)) my_ncol = ncol
675 IF (dbcsr_get_matrix_type(matrix_a) .EQ. dbcsr_type_symmetric) my_alpha2 = my_alpha
676 IF (dbcsr_get_matrix_type(matrix_a) .EQ. dbcsr_type_antisymmetric) my_alpha2 = -my_alpha
678 has_symm = (dbcsr_get_matrix_type(matrix_a) .EQ. dbcsr_type_symmetric .OR. &
679 dbcsr_get_matrix_type(matrix_a) .EQ. dbcsr_type_antisymmetric)
683 CALL dbcsr_iterator_start(iter, matrix_a, read_only=.true., dynamic=.true., dynamic_byrows=.true.)
684 DO WHILE (dbcsr_iterator_blocks_left(iter))
685 CALL dbcsr_iterator_next_block(iter, row, col, data_d, blk, row_offset=rowoff, col_offset=coloff)
686 IF (my_ncol .NE. 1)
THEN
687 CALL dgemm(
'N',
'N', &
688 SIZE(data_d, 1), my_ncol,
SIZE(data_d, 2), &
689 my_alpha, data_d(1, 1),
SIZE(data_d, 1), &
690 vec_b(coloff, 1),
SIZE(vec_b, 1), &
691 1.0_dp, vec_c(rowoff, 1),
SIZE(vec_c, 1))
693 CALL dgemv(
'N',
SIZE(data_d, 1),
SIZE(data_d, 2), &
694 my_alpha, data_d(1, 1),
SIZE(data_d, 1), &
695 vec_b(coloff, 1), 1, &
696 1.0_dp, vec_c(rowoff, 1), 1)
699 CALL dbcsr_iterator_stop(iter)
705 CALL dbcsr_iterator_start(iter, matrix_a)
706 DO WHILE (dbcsr_iterator_blocks_left(iter))
707 CALL dbcsr_iterator_next_block(iter, row, col, data_d, blk, row_offset=rowoff, col_offset=coloff)
708 IF (row .NE. col)
THEN
709 IF (my_ncol .NE. 1)
THEN
710 CALL dgemm(
'T',
'N', &
711 SIZE(data_d, 2), my_ncol,
SIZE(data_d, 1), &
712 my_alpha2, data_d(1, 1),
SIZE(data_d, 1), &
713 vec_b(rowoff, 1),
SIZE(vec_b, 1), &
714 1.0_dp, vec_c(coloff, 1),
SIZE(vec_c, 1))
716 CALL dgemv(
'T',
SIZE(data_d, 1),
SIZE(data_d, 2), &
717 my_alpha2, data_d(1, 1),
SIZE(data_d, 1), &
718 vec_b(rowoff, 1), 1, &
719 1.0_dp, vec_c(coloff, 1), 1)
723 CALL dbcsr_iterator_stop(iter)
726 CALL timestop(timing_handle)
744 TYPE(dbcsr_type),
INTENT(IN) :: matrix
745 TYPE(cp_fm_type),
INTENT(IN) :: fm_in, fm_out
746 INTEGER,
INTENT(IN) :: ncol
747 REAL(
dp),
INTENT(IN),
OPTIONAL :: alpha, beta
749 CHARACTER(LEN=*),
PARAMETER :: routinen =
'cp_dbcsr_sm_fm_multiply'
751 INTEGER :: k_in, k_out, timing_handle, &
752 timing_handle_mult, &
753 a_ncol, a_nrow, b_ncol, b_nrow, c_ncol, c_nrow
754 INTEGER,
DIMENSION(:),
POINTER :: col_blk_size_right_in, &
755 col_blk_size_right_out, row_blk_size, &
757 row_dist, col_dist, col_blk_size
758 TYPE(dbcsr_type) :: in, out
759 REAL(
dp) :: my_alpha, my_beta
760 TYPE(dbcsr_distribution_type) :: dist, dist_right_in, product_dist
762 CALL timeset(routinen, timing_handle)
766 IF (
PRESENT(alpha)) my_alpha = alpha
767 IF (
PRESENT(beta)) my_beta = beta
770 CALL cp_fm_get_info(fm_in, ncol_global=b_ncol, nrow_global=b_nrow)
771 CALL cp_fm_get_info(fm_out, ncol_global=c_ncol, nrow_global=c_nrow)
772 CALL dbcsr_get_info(matrix, nfullrows_total=a_nrow, nfullcols_total=a_ncol)
786 IF (ncol .GT. 0 .AND. k_out .GT. 0 .AND. k_in .GT. 0)
THEN
787 CALL dbcsr_get_info(matrix, row_blk_size=row_blk_size, col_blk_size=col_blk_size, distribution=dist)
790 CALL dbcsr_create(in,
"D", dist_right_in, dbcsr_type_no_symmetry, &
791 col_blk_size, col_blk_size_right_in, nze=0)
793 CALL dbcsr_distribution_get(dist, row_dist=row_dist)
794 CALL dbcsr_distribution_get(dist_right_in, col_dist=col_dist)
795 CALL dbcsr_distribution_new(product_dist, template=dist, &
796 row_dist=row_dist, col_dist=col_dist)
797 ALLOCATE (col_blk_size_right_out(
SIZE(col_blk_size_right_in)))
798 col_blk_size_right_out = col_blk_size_right_in
799 CALL match_col_sizes(col_blk_size_right_out, col_blk_size_right_in, k_out)
806 CALL dbcsr_create(out,
"D", product_dist, dbcsr_type_no_symmetry, &
807 row_blk_size, col_blk_size_right_out, nze=0)
810 IF (ncol .NE. k_out .OR. my_beta .NE. 0.0_dp) &
813 CALL timeset(routinen//
'_core', timing_handle_mult)
814 CALL dbcsr_multiply(
"N",
"N", my_alpha, matrix, in, my_beta, out, &
816 CALL timestop(timing_handle_mult)
820 CALL dbcsr_release(in)
821 CALL dbcsr_release(out)
822 DEALLOCATE (col_blk_size_right_in, col_blk_size_right_out)
823 CALL dbcsr_distribution_release(dist_right_in)
824 CALL dbcsr_distribution_release(product_dist)
828 CALL timestop(timing_handle)
838 SUBROUTINE match_col_sizes(sizes1, sizes2, full_num)
839 INTEGER,
DIMENSION(:),
INTENT(INOUT) :: sizes1
840 INTEGER,
DIMENSION(:),
INTENT(IN) :: sizes2
841 INTEGER,
INTENT(IN) :: full_num
843 INTEGER :: left, n1, n2, p, rm, used
848 cpabort(
"distributions must be equal!")
849 sizes1(1:n1) = sizes2(1:n1)
850 used = sum(sizes1(1:n1))
854 IF (used .LT. full_num)
THEN
855 sizes1(n1) = sizes1(n1) + full_num - used
857 left = used - full_num
859 DO WHILE (left .GT. 0 .AND. p .GT. 0)
860 rm = min(left, sizes1(p))
861 sizes1(p) = sizes1(p) - rm
866 END SUBROUTINE match_col_sizes
887 TYPE(dbcsr_type),
INTENT(INOUT) :: sparse_matrix
888 TYPE(cp_fm_type),
INTENT(IN) :: matrix_v
889 TYPE(cp_fm_type),
OPTIONAL,
INTENT(IN) :: matrix_g
890 INTEGER,
INTENT(IN) :: ncol
891 REAL(kind=
dp),
INTENT(IN),
OPTIONAL :: alpha
892 LOGICAL,
INTENT(IN),
OPTIONAL :: keep_sparsity
893 INTEGER,
INTENT(IN),
OPTIONAL :: symmetry_mode
895 CHARACTER(LEN=*),
PARAMETER :: routinen =
'cp_dbcsr_plus_fm_fm_t_native'
897 INTEGER :: npcols, k, nao, timing_handle, data_type, my_symmetry_mode
898 INTEGER,
DIMENSION(:),
POINTER :: col_blk_size_left, &
899 col_dist_left, row_blk_size, row_dist
900 LOGICAL :: check_product, my_keep_sparsity
901 REAL(kind=
dp) :: my_alpha, norm
902 TYPE(dbcsr_type) :: mat_g, mat_v, sparse_matrix2, &
904 TYPE(cp_fm_struct_type),
POINTER :: fm_struct_tmp
905 TYPE(cp_fm_type) :: fm_matrix
906 TYPE(dbcsr_distribution_type) :: dist_left, sparse_dist
908 check_product = .false.
910 CALL timeset(routinen, timing_handle)
912 my_keep_sparsity = .true.
913 IF (
PRESENT(keep_sparsity)) my_keep_sparsity = keep_sparsity
916 IF (
PRESENT(symmetry_mode)) my_symmetry_mode = symmetry_mode
918 NULLIFY (col_dist_left)
920 IF (ncol .GT. 0)
THEN
921 IF (.NOT. dbcsr_valid_index(sparse_matrix)) &
922 cpabort(
"sparse_matrix must pre-exist")
927 CALL dbcsr_get_info(sparse_matrix, distribution=sparse_dist)
928 CALL dbcsr_distribution_get(sparse_dist, npcols=npcols, row_dist=row_dist)
929 CALL create_bl_distribution(col_dist_left, col_blk_size_left, k, npcols)
930 CALL dbcsr_distribution_new(dist_left, template=sparse_dist, &
931 row_dist=row_dist, col_dist=col_dist_left)
932 DEALLOCATE (col_dist_left)
933 CALL dbcsr_get_info(sparse_matrix, row_blk_size=row_blk_size, data_type=data_type)
934 CALL dbcsr_create(mat_v,
"DBCSR matrix_v", dist_left, dbcsr_type_no_symmetry, &
935 row_blk_size, col_blk_size_left, nze=0, data_type=data_type)
937 CALL dbcsr_verify_matrix(mat_v)
940 IF (
PRESENT(matrix_g))
THEN
941 CALL dbcsr_create(mat_g,
"DBCSR matrix_g", dist_left, dbcsr_type_no_symmetry, &
942 row_blk_size, col_blk_size_left, data_type=data_type)
946 DEALLOCATE (col_blk_size_left)
947 CALL dbcsr_distribution_release(dist_left)
950 IF (check_product)
THEN
952 CALL cp_fm_struct_create(fm_struct_tmp, context=matrix_v%matrix_struct%context, nrow_global=nao, &
953 ncol_global=nao, para_env=matrix_v%matrix_struct%para_env)
954 CALL cp_fm_create(fm_matrix, fm_struct_tmp, name=
"fm matrix")
957 CALL dbcsr_copy(sparse_matrix3, sparse_matrix)
961 IF (
PRESENT(alpha)) my_alpha = alpha
962 IF (
PRESENT(matrix_g))
THEN
963 IF (my_symmetry_mode == 1)
THEN
965 CALL dbcsr_multiply(
"N",
"T", 0.5_dp*my_alpha, mat_v, mat_g, &
966 1.0_dp, sparse_matrix, &
967 retain_sparsity=my_keep_sparsity, &
969 CALL dbcsr_multiply(
"N",
"T", 0.5_dp*my_alpha, mat_g, mat_v, &
970 1.0_dp, sparse_matrix, &
971 retain_sparsity=my_keep_sparsity, &
973 ELSE IF (my_symmetry_mode == -1)
THEN
975 CALL dbcsr_multiply(
"N",
"T", 0.5_dp*my_alpha, mat_v, mat_g, &
976 1.0_dp, sparse_matrix, &
977 retain_sparsity=my_keep_sparsity, &
979 CALL dbcsr_multiply(
"N",
"T", -0.5_dp*my_alpha, mat_g, mat_v, &
980 1.0_dp, sparse_matrix, &
981 retain_sparsity=my_keep_sparsity, &
985 CALL dbcsr_multiply(
"N",
"T", my_alpha, mat_v, mat_g, &
986 1.0_dp, sparse_matrix, &
987 retain_sparsity=my_keep_sparsity, &
991 CALL dbcsr_multiply(
"N",
"T", my_alpha, mat_v, mat_v, &
992 1.0_dp, sparse_matrix, &
993 retain_sparsity=my_keep_sparsity, &
997 IF (check_product)
THEN
998 IF (
PRESENT(matrix_g))
THEN
999 IF (my_symmetry_mode == 1)
THEN
1000 CALL cp_fm_gemm(
"N",
"T", nao, nao, ncol, 0.5_dp*my_alpha, matrix_v, matrix_g, &
1002 CALL cp_fm_gemm(
"N",
"T", nao, nao, ncol, 0.5_dp*my_alpha, matrix_g, matrix_v, &
1004 ELSE IF (my_symmetry_mode == -1)
THEN
1005 CALL cp_fm_gemm(
"N",
"T", nao, nao, ncol, 0.5_dp*my_alpha, matrix_v, matrix_g, &
1007 CALL cp_fm_gemm(
"N",
"T", nao, nao, ncol, -0.5_dp*my_alpha, matrix_g, matrix_v, &
1010 CALL cp_fm_gemm(
"N",
"T", nao, nao, ncol, my_alpha, matrix_v, matrix_g, &
1014 CALL cp_fm_gemm(
"N",
"T", nao, nao, ncol, my_alpha, matrix_v, matrix_v, &
1018 CALL dbcsr_copy(sparse_matrix2, sparse_matrix)
1019 CALL dbcsr_scale(sparse_matrix2, alpha_scalar=0.0_dp)
1020 CALL copy_fm_to_dbcsr(fm_matrix, sparse_matrix2, keep_sparsity=my_keep_sparsity)
1021 CALL dbcsr_add(sparse_matrix2, sparse_matrix, alpha_scalar=1.0_dp, &
1022 beta_scalar=-1.0_dp)
1023 CALL dbcsr_norm(sparse_matrix2, which_norm=dbcsr_norm_frobenius, &
1025 WRITE (*, *)
'nao=', nao,
' k=', k,
' ncol=', ncol,
' my_alpha=', my_alpha
1026 WRITE (*, *)
'PRESENT (matrix_g)',
PRESENT(matrix_g)
1027 WRITE (*, *)
'matrix_type=', dbcsr_get_matrix_type(sparse_matrix)
1028 WRITE (*, *)
'norm(sm+alpha*v*g^t - fm+alpha*v*g^t)/n=', norm/real(nao,
dp)
1029 IF (norm/real(nao,
dp) .GT. 1e-12_dp)
THEN
1048 CALL dbcsr_release(sparse_matrix2)
1049 CALL dbcsr_release(sparse_matrix3)
1050 CALL cp_fm_release(fm_matrix)
1052 CALL dbcsr_release(mat_v)
1053 IF (
PRESENT(matrix_g))
CALL dbcsr_release(mat_g)
1055 CALL timestop(timing_handle)
1070 TYPE(dbcsr_type),
INTENT(INOUT) :: matrix
1071 TYPE(cp_fm_type),
INTENT(IN) :: fm_in
1072 TYPE(dbcsr_type),
INTENT(IN) :: template
1074 INTEGER :: k_in, data_type
1075 INTEGER,
DIMENSION(:),
POINTER :: col_blk_size_right_in, row_blk_size
1076 TYPE(dbcsr_distribution_type) :: tmpl_dist, dist_right_in
1080 CALL dbcsr_get_info(template, distribution=tmpl_dist)
1082 CALL dbcsr_get_info(template, row_blk_size=row_blk_size, data_type=data_type)
1083 CALL dbcsr_create(matrix,
"D", dist_right_in, dbcsr_type_no_symmetry, &
1084 row_blk_size, col_blk_size_right_in, nze=0, data_type=data_type)
1087 DEALLOCATE (col_blk_size_right_in)
1088 CALL dbcsr_distribution_release(dist_right_in)
1105 TYPE(dbcsr_type),
INTENT(INOUT) :: matrix, template
1106 INTEGER,
INTENT(IN) :: m, n
1107 CHARACTER,
OPTIONAL,
INTENT(IN) :: sym
1108 INTEGER,
OPTIONAL,
INTENT(IN) :: data_type
1111 INTEGER :: my_data_type, nprows, npcols
1112 INTEGER,
DIMENSION(:),
POINTER :: col_blk_size, &
1113 col_dist, row_blk_size, &
1115 TYPE(dbcsr_distribution_type) :: tmpl_dist, dist_m_n
1117 CALL dbcsr_get_info(template, &
1118 matrix_type=mysym, &
1119 data_type=my_data_type, &
1120 distribution=tmpl_dist)
1122 IF (
PRESENT(sym)) mysym = sym
1123 IF (
PRESENT(data_type)) my_data_type = data_type
1125 NULLIFY (row_dist, col_dist)
1126 NULLIFY (row_blk_size, col_blk_size)
1129 CALL dbcsr_distribution_get(tmpl_dist, nprows=nprows, npcols=npcols)
1130 CALL create_bl_distribution(row_dist, row_blk_size, m, nprows)
1131 CALL create_bl_distribution(col_dist, col_blk_size, n, npcols)
1132 CALL dbcsr_distribution_new(dist_m_n, template=tmpl_dist, &
1133 row_dist=row_dist, col_dist=col_dist, &
1135 reuse_arrays=.true.)
1137 CALL dbcsr_create(matrix,
"m_n_template", dist_m_n, mysym, &
1138 row_blk_size, col_blk_size, nze=0, data_type=my_data_type, &
1139 reuse_arrays=.true.)
1140 CALL dbcsr_distribution_release(dist_m_n)
1156 TYPE(dbcsr_type),
INTENT(INOUT) :: matrix, template
1158 CHARACTER,
OPTIONAL :: sym
1159 INTEGER,
OPTIONAL :: data_type
1162 INTEGER :: my_data_type, npcols
1163 INTEGER,
DIMENSION(:),
POINTER :: col_blk_size, col_dist, row_blk_size, &
1165 TYPE(dbcsr_distribution_type) :: dist_m_n, tmpl_dist
1167 mysym = dbcsr_get_matrix_type(template)
1168 IF (
PRESENT(sym)) mysym = sym
1169 my_data_type = dbcsr_get_data_type(template)
1170 IF (
PRESENT(data_type)) my_data_type = data_type
1172 CALL dbcsr_get_info(template, distribution=tmpl_dist)
1173 CALL dbcsr_distribution_get(tmpl_dist, &
1177 NULLIFY (col_dist, col_blk_size)
1178 CALL create_bl_distribution(col_dist, col_blk_size, n, npcols)
1179 CALL dbcsr_distribution_new(dist_m_n, template=tmpl_dist, &
1180 row_dist=row_dist, col_dist=col_dist)
1182 CALL dbcsr_get_info(template, row_blk_size=row_blk_size)
1183 CALL dbcsr_create(matrix,
"m_n_template", dist_m_n, mysym, &
1184 row_blk_size, col_blk_size, nze=0, data_type=my_data_type)
1186 DEALLOCATE (col_dist, col_blk_size)
1187 CALL dbcsr_distribution_release(dist_m_n)
1201 SUBROUTINE create_bl_distribution(block_distribution, &
1202 block_size, nelements, nbins)
1203 INTEGER,
DIMENSION(:),
INTENT(OUT),
POINTER :: block_distribution, block_size
1204 INTEGER,
INTENT(IN) :: nelements, nbins
1206 CHARACTER(len=*),
PARAMETER :: routinen =
'create_bl_distribution', &
1207 routinep = modulen//
':'//routinen
1209 INTEGER :: bin, blk_layer, element_stack, els, &
1210 estimated_blocks, max_blocks_per_bin, &
1211 nblks, nblocks, stat
1212 INTEGER,
DIMENSION(:),
POINTER :: blk_dist, blk_sizes
1216 NULLIFY (block_distribution)
1217 NULLIFY (block_size)
1219 IF (nelements .GT. 0)
THEN
1222 max_blocks_per_bin = ceiling(real(nblocks, kind=
dp)/real(nbins, kind=
dp))
1225 WRITE (*,
'(1X,A,1X,A,I7,A,I7,A)') routinep,
"For", nelements, &
1226 " elements and", nbins,
" bins"
1227 WRITE (*,
'(1X,A,1X,A,I7,A)') routinep,
"There are", &
1229 WRITE (*,
'(1X,A,1X,A,I7,A)') routinep,
"There are", &
1231 WRITE (*,
'(1X,A,1X,A,I7,A)') routinep,
"There are", &
1232 max_blocks_per_bin,
" max blocks/bin"
1235 estimated_blocks = max_blocks_per_bin*nbins
1236 ALLOCATE (blk_dist(estimated_blocks), stat=stat)
1239 ALLOCATE (blk_sizes(estimated_blocks), stat=stat)
1241 cpabort(
"blk_sizes")
1244 DO blk_layer = 1, max_blocks_per_bin
1245 DO bin = 0, nbins - 1
1247 IF (els .GT. 0)
THEN
1248 element_stack = element_stack + els
1250 blk_dist(nblks) = bin
1251 blk_sizes(nblks) = els
1252 IF (
debug_mod)
WRITE (*,
'(1X,A,I5,A,I5,A,I5)') routinep//
" Assigning", &
1253 els,
" elements as block", nblks,
" to bin", bin
1258 IF (nblks .EQ. estimated_blocks)
THEN
1259 block_distribution => blk_dist
1260 block_size => blk_sizes
1262 ALLOCATE (block_distribution(nblks), stat=stat)
1265 block_distribution(:) = blk_dist(1:nblks)
1266 DEALLOCATE (blk_dist)
1267 ALLOCATE (block_size(nblks), stat=stat)
1269 cpabort(
"blk_sizes")
1270 block_size(:) = blk_sizes(1:nblks)
1271 DEALLOCATE (blk_sizes)
1274 ALLOCATE (block_distribution(0), stat=stat)
1277 ALLOCATE (block_size(0), stat=stat)
1279 cpabort(
"blk_sizes")
1281 1579
FORMAT(i5, 1x, i5, 1x, i5, 1x, i5, 1x, i5, 1x, i5, 1x, i5, 1x, i5, 1x, i5, 1x, i5)
1283 WRITE (*,
'(1X,A,A)') routinep//
" Distribution"
1284 WRITE (*, 1579) block_distribution(:)
1285 WRITE (*,
'(1X,A,A)') routinep//
" Sizes"
1286 WRITE (*, 1579) block_size(:)
1288 END SUBROUTINE create_bl_distribution
1302 right_col_blk_sizes)
1303 TYPE(dbcsr_distribution_type),
INTENT(OUT) :: dist_right
1304 TYPE(dbcsr_distribution_type),
INTENT(IN) :: dist_left
1305 INTEGER,
INTENT(IN) :: ncolumns
1306 INTEGER,
DIMENSION(:),
INTENT(OUT),
POINTER :: right_col_blk_sizes
1308 INTEGER :: multiplicity, nimages, ncols, nprows, npcols
1309 INTEGER,
ALLOCATABLE,
DIMENSION(:) :: tmp_images
1310 INTEGER,
DIMENSION(:),
POINTER :: right_col_dist, right_row_dist, old_col_dist, dummy
1312 CALL dbcsr_distribution_get(dist_left, &
1314 col_dist=old_col_dist, &
1319 CALL create_bl_distribution(right_col_dist, right_col_blk_sizes, ncolumns, npcols)
1321 ALLOCATE (right_row_dist(ncols), tmp_images(ncols))
1322 nimages =
lcm(nprows, npcols)/nprows
1323 multiplicity = nprows/
gcd(nprows, npcols)
1324 CALL rebin_distribution(right_row_dist, tmp_images, old_col_dist, nprows, multiplicity, nimages)
1327 CALL dbcsr_distribution_new(dist_right, &
1328 template=dist_left, &
1329 row_dist=right_row_dist, &
1330 col_dist=right_col_dist, &
1333 reuse_arrays=.true.)
1334 DEALLOCATE (tmp_images)
1370 SUBROUTINE rebin_distribution(new_bins, images, source_bins, &
1371 nbins, multiplicity, nimages)
1372 INTEGER,
DIMENSION(:),
INTENT(OUT) :: new_bins, images
1373 INTEGER,
DIMENSION(:),
INTENT(IN) :: source_bins
1374 INTEGER,
INTENT(IN) :: nbins, multiplicity, nimages
1376 INTEGER :: bin, i, old_nbins, virtual_bin
1377 INTEGER,
ALLOCATABLE,
DIMENSION(:) :: bin_multiplier
1381 IF (mod(nbins*nimages, multiplicity) .NE. 0) &
1382 cpwarn(
"mulitplicity is not divisor of new process grid coordinate")
1383 old_nbins = (nbins*nimages)/multiplicity
1384 ALLOCATE (bin_multiplier(0:old_nbins - 1))
1385 bin_multiplier(:) = 0
1386 DO i = 1,
SIZE(new_bins)
1387 IF (i .LE.
SIZE(source_bins))
THEN
1388 bin = source_bins(i)
1391 bin = mod(i, old_nbins)
1393 virtual_bin = bin*multiplicity + bin_multiplier(bin)
1394 new_bins(i) = virtual_bin/nimages
1395 images(i) = 1 + mod(virtual_bin, nimages)
1396 bin_multiplier(bin) = bin_multiplier(bin) + 1
1397 IF (bin_multiplier(bin) .GE. multiplicity)
THEN
1398 bin_multiplier(bin) = 0
1401 END SUBROUTINE rebin_distribution
1417 SUBROUTINE dbcsr_create_dist_block_cyclic(dist, nrows, ncolumns, &
1418 nrow_block, ncol_block, group_handle, pgrid, row_blk_sizes, col_blk_sizes)
1419 TYPE(dbcsr_distribution_type),
INTENT(OUT) :: dist
1420 INTEGER,
INTENT(IN) :: nrows, ncolumns, nrow_block, ncol_block
1421 INTEGER,
INTENT(IN) :: group_handle
1422 INTEGER,
DIMENSION(:, :),
POINTER :: pgrid
1423 INTEGER,
DIMENSION(:),
INTENT(OUT),
POINTER :: row_blk_sizes, col_blk_sizes
1425 CHARACTER(len=*),
PARAMETER :: routinen =
'dbcsr_create_dist_block_cyclic'
1427 INTEGER :: nblkcols, nblkrows, npcols, nprows, &
1429 INTEGER,
DIMENSION(:),
POINTER :: cd_data, rd_data
1432 IF (nrow_block .EQ. 0)
THEN
1436 nblkrows = nrows/nrow_block
1437 sz = mod(nrows, nrow_block)
1439 IF (sz .GT. 0) nblkrows = nblkrows + 1
1440 ALLOCATE (row_blk_sizes(nblkrows), rd_data(nblkrows))
1441 row_blk_sizes = nrow_block
1442 IF (sz .NE. 0) row_blk_sizes(nblkrows) = sz
1445 IF (ncol_block .EQ. 0)
THEN
1449 nblkcols = ncolumns/ncol_block
1450 sz = mod(ncolumns, ncol_block)
1452 IF (sz .GT. 0) nblkcols = nblkcols + 1
1453 ALLOCATE (col_blk_sizes(nblkcols), cd_data(nblkcols))
1454 col_blk_sizes = ncol_block
1455 IF (sz .NE. 0) col_blk_sizes(nblkcols) = sz
1458 WRITE (*, *) routinen//
" nrows,nrow_block,nblkrows=", &
1459 nrows, nrow_block, nblkrows
1460 WRITE (*, *) routinen//
" ncols,ncol_block,nblkcols=", &
1461 ncolumns, ncol_block, nblkcols
1464 nprows =
SIZE(pgrid, 1)
1465 DO pdim = 0, min(nprows - 1, nblkrows - 1)
1466 rd_data(1 + pdim:nblkrows:nprows) = pdim
1469 npcols =
SIZE(pgrid, 2)
1470 DO pdim = 0, min(npcols - 1, nblkcols - 1)
1471 cd_data(1 + pdim:nblkcols:npcols) = pdim
1475 WRITE (*, *) routinen//
" row_dist", &
1477 WRITE (*, *) routinen//
" col_dist", &
1481 CALL dbcsr_distribution_new(dist, &
1482 group=group_handle, pgrid=pgrid, &
1485 reuse_arrays=.true.)
1487 END SUBROUTINE dbcsr_create_dist_block_cyclic
1496 SUBROUTINE allocate_dbcsr_matrix_set_1d(matrix_set, nmatrix)
1497 TYPE(dbcsr_p_type),
DIMENSION(:),
POINTER :: matrix_set
1498 INTEGER,
INTENT(IN) :: nmatrix
1503 ALLOCATE (matrix_set(nmatrix))
1504 DO imatrix = 1, nmatrix
1505 NULLIFY (matrix_set(imatrix)%matrix)
1507 END SUBROUTINE allocate_dbcsr_matrix_set_1d
1517 SUBROUTINE allocate_dbcsr_matrix_set_2d(matrix_set, nmatrix, mmatrix)
1518 TYPE(dbcsr_p_type),
DIMENSION(:, :),
POINTER :: matrix_set
1519 INTEGER,
INTENT(IN) :: nmatrix, mmatrix
1521 INTEGER :: imatrix, jmatrix
1524 ALLOCATE (matrix_set(nmatrix, mmatrix))
1525 DO jmatrix = 1, mmatrix
1526 DO imatrix = 1, nmatrix
1527 NULLIFY (matrix_set(imatrix, jmatrix)%matrix)
1530 END SUBROUTINE allocate_dbcsr_matrix_set_2d
1541 SUBROUTINE allocate_dbcsr_matrix_set_3d(matrix_set, nmatrix, mmatrix, pmatrix)
1542 TYPE(dbcsr_p_type),
DIMENSION(:, :, :),
POINTER :: matrix_set
1543 INTEGER,
INTENT(IN) :: nmatrix, mmatrix, pmatrix
1545 INTEGER :: imatrix, jmatrix, kmatrix
1548 ALLOCATE (matrix_set(nmatrix, mmatrix, pmatrix))
1549 DO kmatrix = 1, pmatrix
1550 DO jmatrix = 1, mmatrix
1551 DO imatrix = 1, nmatrix
1552 NULLIFY (matrix_set(imatrix, jmatrix, kmatrix)%matrix)
1556 END SUBROUTINE allocate_dbcsr_matrix_set_3d
1567 SUBROUTINE allocate_dbcsr_matrix_set_4d(matrix_set, nmatrix, mmatrix, pmatrix, qmatrix)
1568 TYPE(dbcsr_p_type),
DIMENSION(:, :, :, :),
POINTER :: matrix_set
1569 INTEGER,
INTENT(IN) :: nmatrix, mmatrix, pmatrix, qmatrix
1571 INTEGER :: imatrix, jmatrix, kmatrix, lmatrix
1574 ALLOCATE (matrix_set(nmatrix, mmatrix, pmatrix, qmatrix))
1575 DO lmatrix = 1, qmatrix
1576 DO kmatrix = 1, pmatrix
1577 DO jmatrix = 1, mmatrix
1578 DO imatrix = 1, nmatrix
1579 NULLIFY (matrix_set(imatrix, jmatrix, kmatrix, lmatrix)%matrix)
1584 END SUBROUTINE allocate_dbcsr_matrix_set_4d
1595 SUBROUTINE allocate_dbcsr_matrix_set_5d(matrix_set, nmatrix, mmatrix, pmatrix, qmatrix, smatrix)
1596 TYPE(dbcsr_p_type),
DIMENSION(:, :, :, :, :),
POINTER :: matrix_set
1597 INTEGER,
INTENT(IN) :: nmatrix, mmatrix, pmatrix, qmatrix, smatrix
1599 INTEGER :: imatrix, jmatrix, kmatrix, lmatrix, hmatrix
1602 ALLOCATE (matrix_set(nmatrix, mmatrix, pmatrix, qmatrix, smatrix))
1603 DO hmatrix = 1, smatrix
1604 DO lmatrix = 1, qmatrix
1605 DO kmatrix = 1, pmatrix
1606 DO jmatrix = 1, mmatrix
1607 DO imatrix = 1, nmatrix
1608 NULLIFY (matrix_set(imatrix, jmatrix, kmatrix, lmatrix, hmatrix)%matrix)
1614 END SUBROUTINE allocate_dbcsr_matrix_set_5d
1622 SUBROUTINE deallocate_dbcsr_matrix_set_1d(matrix_set)
1624 TYPE(dbcsr_p_type),
DIMENSION(:),
POINTER :: matrix_set
1628 IF (
ASSOCIATED(matrix_set))
THEN
1629 DO imatrix = 1,
SIZE(matrix_set)
1630 CALL dbcsr_deallocate_matrix(matrix_set(imatrix)%matrix)
1632 DEALLOCATE (matrix_set)
1635 END SUBROUTINE deallocate_dbcsr_matrix_set_1d
1643 SUBROUTINE deallocate_dbcsr_matrix_set_2d(matrix_set)
1645 TYPE(dbcsr_p_type),
DIMENSION(:, :),
POINTER :: matrix_set
1647 INTEGER :: imatrix, jmatrix
1649 IF (
ASSOCIATED(matrix_set))
THEN
1650 DO jmatrix = 1,
SIZE(matrix_set, 2)
1651 DO imatrix = 1,
SIZE(matrix_set, 1)
1652 CALL dbcsr_deallocate_matrix(matrix_set(imatrix, jmatrix)%matrix)
1655 DEALLOCATE (matrix_set)
1657 END SUBROUTINE deallocate_dbcsr_matrix_set_2d
1665 SUBROUTINE deallocate_dbcsr_matrix_set_3d(matrix_set)
1667 TYPE(dbcsr_p_type),
DIMENSION(:, :, :),
POINTER :: matrix_set
1669 INTEGER :: imatrix, jmatrix, kmatrix
1671 IF (
ASSOCIATED(matrix_set))
THEN
1672 DO kmatrix = 1,
SIZE(matrix_set, 3)
1673 DO jmatrix = 1,
SIZE(matrix_set, 2)
1674 DO imatrix = 1,
SIZE(matrix_set, 1)
1675 CALL dbcsr_deallocate_matrix(matrix_set(imatrix, jmatrix, kmatrix)%matrix)
1679 DEALLOCATE (matrix_set)
1681 END SUBROUTINE deallocate_dbcsr_matrix_set_3d
1689 SUBROUTINE deallocate_dbcsr_matrix_set_4d(matrix_set)
1691 TYPE(dbcsr_p_type),
DIMENSION(:, :, :, :),
POINTER :: matrix_set
1693 INTEGER :: imatrix, jmatrix, kmatrix, lmatrix
1695 IF (
ASSOCIATED(matrix_set))
THEN
1696 DO lmatrix = 1,
SIZE(matrix_set, 4)
1697 DO kmatrix = 1,
SIZE(matrix_set, 3)
1698 DO jmatrix = 1,
SIZE(matrix_set, 2)
1699 DO imatrix = 1,
SIZE(matrix_set, 1)
1700 CALL dbcsr_deallocate_matrix(matrix_set(imatrix, jmatrix, kmatrix, lmatrix)%matrix)
1705 DEALLOCATE (matrix_set)
1707 END SUBROUTINE deallocate_dbcsr_matrix_set_4d
1715 SUBROUTINE deallocate_dbcsr_matrix_set_5d(matrix_set)
1717 TYPE(dbcsr_p_type),
DIMENSION(:, :, :, :, :),
POINTER :: matrix_set
1719 INTEGER :: imatrix, jmatrix, kmatrix, hmatrix, lmatrix
1721 IF (
ASSOCIATED(matrix_set))
THEN
1722 DO hmatrix = 1,
SIZE(matrix_set, 5)
1723 DO lmatrix = 1,
SIZE(matrix_set, 4)
1724 DO kmatrix = 1,
SIZE(matrix_set, 3)
1725 DO jmatrix = 1,
SIZE(matrix_set, 2)
1726 DO imatrix = 1,
SIZE(matrix_set, 1)
1727 CALL dbcsr_deallocate_matrix(matrix_set(imatrix, jmatrix, kmatrix, lmatrix, hmatrix)%matrix)
1733 DEALLOCATE (matrix_set)
1735 END SUBROUTINE deallocate_dbcsr_matrix_set_5d
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.
methods related to the blacs parallel environment
Represents a complex full matrix distributed on many processors.
DBCSR operations in CP2K.
subroutine, public cp_dbcsr_sm_fm_multiply(matrix, fm_in, fm_out, ncol, alpha, beta)
multiply a dbcsr with a fm matrix
integer, save, public max_elements_per_block
subroutine, public copy_cfm_to_dbcsr(fm, matrix, keep_sparsity)
Copy a BLACS matrix to a dbcsr matrix.
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.
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
logical, parameter debug_mod
subroutine, public cp_dbcsr_dist2d_to_dist(dist2d, dist)
Creates a DBCSR distribution from a distribution_2d.
subroutine, public cp_dbcsr_m_by_n_from_template(matrix, template, m, n, sym, data_type)
Utility function to create an arbitrary shaped dbcsr matrix with the same processor grid as the templ...
subroutine, public copy_dbcsr_to_fm(matrix, fm)
Copy a DBCSR matrix to a BLACS matrix.
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...
subroutine, public copy_dbcsr_to_fm_bc(bc_mat, fm)
Copy a DBCSR_BLACS matrix to a BLACS matrix.
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...
subroutine, public copy_dbcsr_to_cfm(matrix, fm)
Copy a DBCSR matrix to a BLACS matrix.
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,...
subroutine, public cp_dbcsr_m_by_n_from_row_template(matrix, template, n, sym, data_type)
Utility function to create dbcsr matrix, m x n matrix (n arbitrary) with the same processor grid and ...
subroutine, public dbcsr_copy_columns_hack(matrix_b, matrix_a, ncol, source_start, target_start, para_env, blacs_env)
hack for dbcsr_copy_columns
subroutine, public copy_fm_to_dbcsr(fm, matrix, keep_sparsity)
Copy a BLACS matrix to a dbcsr matrix.
basic linear algebra operations for full matrices
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 )
represent the structure of a full matrix
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
subroutine, public cp_fm_struct_release(fmstruct)
releases a full matrix structure
represent a full matrix distributed on many processors
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
subroutine, public cp_fm_create(matrix, matrix_struct, name, use_sp)
creates a new full matrix with the given structure
stores a mapping of 2D info (e.g. matrix) on a 2D processor distribution (i.e. blacs grid) where cpus...
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
Defines the basic variable types.
integer, parameter, public dp
integer, parameter, public default_string_length
Collection of simple mathematical functions and subroutines.
elemental integer function, public lcm(a, b)
computes the least common multiplier of two numbers
elemental integer function, public gcd(a, b)
computes the greatest common divisor of two number
Interface to the message passing library MPI.