19 USE dbcsr_api,
ONLY: &
20 dbcsr_distribution_get, dbcsr_distribution_type, dbcsr_filter, dbcsr_finalize, &
21 dbcsr_get_block_p, dbcsr_get_info, dbcsr_get_matrix_type, dbcsr_get_stored_coordinates, &
22 dbcsr_iterator_blocks_left, dbcsr_iterator_next_block, dbcsr_iterator_start, &
23 dbcsr_iterator_stop, dbcsr_iterator_type, dbcsr_nblkcols_total, dbcsr_nblkrows_total, &
24 dbcsr_reserve_block2d, dbcsr_type, dbcsr_type_antisymmetric, dbcsr_type_no_symmetry, &
25 dbcsr_type_symmetric, dbcsr_work_create
27 domain_submatrix_type,&
32 #include "./base/base_uses.f90"
38 CHARACTER(len=*),
PARAMETER,
PRIVATE :: moduleN =
'domain_submatrix_methods'
41 release_submatrices, multiply_submatrices, add_submatrices, &
47 INTERFACE init_submatrices
48 MODULE PROCEDURE init_submatrices_0d
49 MODULE PROCEDURE init_submatrices_1d
50 MODULE PROCEDURE init_submatrices_2d
53 INTERFACE set_submatrices
54 MODULE PROCEDURE set_submatrix_array
55 MODULE PROCEDURE set_submatrix
58 INTERFACE copy_submatrices
59 MODULE PROCEDURE copy_submatrix_array
60 MODULE PROCEDURE copy_submatrix
63 INTERFACE release_submatrices
64 MODULE PROCEDURE release_submatrix_array
65 MODULE PROCEDURE release_submatrix
68 INTERFACE multiply_submatrices
69 MODULE PROCEDURE multiply_submatrices_once
70 MODULE PROCEDURE multiply_submatrices_array
73 INTERFACE add_submatrices
74 MODULE PROCEDURE add_submatrices_once
75 MODULE PROCEDURE add_submatrices_array
84 SUBROUTINE init_submatrices_0d(subm)
86 TYPE(domain_submatrix_type),
INTENT(INOUT) :: subm
96 END SUBROUTINE init_submatrices_0d
102 SUBROUTINE init_submatrices_1d(subm)
104 TYPE(domain_submatrix_type),
DIMENSION(:), &
105 INTENT(INOUT) :: subm
115 END SUBROUTINE init_submatrices_1d
121 SUBROUTINE init_submatrices_2d(subm)
123 TYPE(domain_submatrix_type),
DIMENSION(:, :), &
124 INTENT(INOUT) :: subm
126 subm(:, :)%domain = -1
127 subm(:, :)%nbrows = -1
128 subm(:, :)%nbcols = -1
129 subm(:, :)%nrows = -1
130 subm(:, :)%ncols = -1
131 subm(:, :)%nnodes = -1
134 END SUBROUTINE init_submatrices_2d
142 SUBROUTINE copy_submatrix_array(original, copy, copy_data)
144 TYPE(domain_submatrix_type),
DIMENSION(:), &
145 INTENT(IN) :: original
146 TYPE(domain_submatrix_type),
DIMENSION(:), &
147 INTENT(INOUT) :: copy
148 LOGICAL,
INTENT(IN) :: copy_data
150 CHARACTER(len=*),
PARAMETER :: routineN =
'copy_submatrix_array'
152 INTEGER :: handle, idomain, ndomains, ndomainsB
154 CALL timeset(routinen, handle)
156 ndomains =
SIZE(original)
157 ndomainsb =
SIZE(copy)
158 cpassert(ndomains .EQ. ndomainsb)
159 copy(:)%nnodes = original(:)%nnodes
160 copy(:)%group = original(:)%group
161 DO idomain = 1, ndomains
162 IF (original(idomain)%domain .GT. 0)
THEN
163 CALL copy_submatrix(original(idomain), copy(idomain), copy_data)
167 CALL timestop(handle)
169 END SUBROUTINE copy_submatrix_array
177 SUBROUTINE copy_submatrix(original, copy, copy_data)
179 TYPE(domain_submatrix_type),
INTENT(IN) :: original
180 TYPE(domain_submatrix_type),
INTENT(INOUT) :: copy
181 LOGICAL,
INTENT(IN) :: copy_data
183 CHARACTER(len=*),
PARAMETER :: routineN =
'copy_submatrix'
185 INTEGER :: handle, icol, irow
187 CALL timeset(routinen, handle)
189 copy%domain = original%domain
190 copy%nnodes = original%nnodes
191 copy%group = original%group
193 IF (original%domain .GT. 0)
THEN
195 copy%nbrows = original%nbrows
196 copy%nbcols = original%nbcols
197 copy%nrows = original%nrows
198 copy%ncols = original%ncols
200 IF (.NOT.
ALLOCATED(copy%dbcsr_row))
THEN
201 ALLOCATE (copy%dbcsr_row(original%nbrows))
203 IF (
SIZE(copy%dbcsr_row) .NE.
SIZE(original%dbcsr_row))
THEN
204 DEALLOCATE (copy%dbcsr_row)
205 ALLOCATE (copy%dbcsr_row(original%nbrows))
208 IF (.NOT.
ALLOCATED(copy%dbcsr_col))
THEN
209 ALLOCATE (copy%dbcsr_col(original%nbcols))
211 IF (
SIZE(copy%dbcsr_col) .NE.
SIZE(original%dbcsr_col))
THEN
212 DEALLOCATE (copy%dbcsr_col)
213 ALLOCATE (copy%dbcsr_col(original%nbcols))
216 IF (.NOT.
ALLOCATED(copy%size_brow))
THEN
217 ALLOCATE (copy%size_brow(original%nbrows))
219 IF (
SIZE(copy%size_brow) .NE.
SIZE(original%size_brow))
THEN
220 DEALLOCATE (copy%size_brow)
221 ALLOCATE (copy%size_brow(original%nbrows))
224 IF (.NOT.
ALLOCATED(copy%size_bcol))
THEN
225 ALLOCATE (copy%size_bcol(original%nbcols))
227 IF (
SIZE(copy%size_bcol) .NE.
SIZE(original%size_bcol))
THEN
228 DEALLOCATE (copy%size_bcol)
229 ALLOCATE (copy%size_bcol(original%nbcols))
233 DO irow = 1, original%nbrows
234 copy%dbcsr_row(irow) = original%dbcsr_row(irow)
235 copy%size_brow(irow) = original%size_brow(irow)
238 DO icol = 1, original%nbcols
239 copy%dbcsr_col(icol) = original%dbcsr_col(icol)
240 copy%size_bcol(icol) = original%size_bcol(icol)
249 CALL timestop(handle)
251 END SUBROUTINE copy_submatrix
260 REAL(kind=
dp),
DIMENSION(:, :),
INTENT(IN) :: array
261 TYPE(domain_submatrix_type),
INTENT(INOUT) :: copy
263 CHARACTER(len=*),
PARAMETER :: routinen =
'copy_submatrix_data'
265 INTEGER :: ds1, ds2, handle, ms1, ms2
267 CALL timeset(routinen, handle)
269 cpassert(copy%domain .GT. 0)
274 IF (.NOT.
ALLOCATED(copy%mdata))
THEN
275 ALLOCATE (copy%mdata(ds1, ds2))
277 ms1 =
SIZE(copy%mdata, 1)
278 ms2 =
SIZE(copy%mdata, 2)
279 IF ((ds1 .NE. ms1) .OR. (ds2 .NE. ms2))
THEN
280 DEALLOCATE (copy%mdata)
281 ALLOCATE (copy%mdata(ds1, ds2))
285 copy%mdata(:, :) = array(:, :)
287 CALL timestop(handle)
296 SUBROUTINE set_submatrix_array(submatrices, scalar)
298 TYPE(domain_submatrix_type),
DIMENSION(:), &
299 INTENT(INOUT) :: submatrices
300 REAL(kind=
dp),
INTENT(IN) :: scalar
302 CHARACTER(len=*),
PARAMETER :: routinen =
'set_submatrix_array'
304 INTEGER :: handle, idomain, ndomains
306 CALL timeset(routinen, handle)
308 ndomains =
SIZE(submatrices)
309 DO idomain = 1, ndomains
310 IF (submatrices(idomain)%domain .GT. 0)
THEN
311 CALL set_submatrix(submatrices(idomain), scalar)
315 CALL timestop(handle)
317 END SUBROUTINE set_submatrix_array
324 SUBROUTINE set_submatrix(submatrix, scalar)
326 TYPE(domain_submatrix_type),
INTENT(INOUT) :: submatrix
327 REAL(kind=
dp),
INTENT(IN) :: scalar
329 CHARACTER(len=*),
PARAMETER :: routinen =
'set_submatrix'
331 INTEGER :: ds1, ds2, handle, ms1, ms2
333 CALL timeset(routinen, handle)
335 cpassert(submatrix%domain .GT. 0)
336 cpassert(submatrix%nrows .GT. 0)
337 cpassert(submatrix%ncols .GT. 0)
339 ds1 = submatrix%nrows
340 ds2 = submatrix%ncols
342 IF (.NOT.
ALLOCATED(submatrix%mdata))
THEN
343 ALLOCATE (submatrix%mdata(ds1, ds2))
345 ms1 =
SIZE(submatrix%mdata, 1)
346 ms2 =
SIZE(submatrix%mdata, 2)
347 IF ((ds1 .NE. ms1) .OR. (ds2 .NE. ms2))
THEN
348 DEALLOCATE (submatrix%mdata)
349 ALLOCATE (submatrix%mdata(ds1, ds2))
353 submatrix%mdata(:, :) = scalar
355 CALL timestop(handle)
357 END SUBROUTINE set_submatrix
363 SUBROUTINE release_submatrix_array(subm)
365 TYPE(domain_submatrix_type),
DIMENSION(:), &
366 INTENT(INOUT) :: subm
368 CHARACTER(len=*),
PARAMETER :: routinen =
'release_submatrix_array'
370 INTEGER :: handle, idomain, ndomains
372 CALL timeset(routinen, handle)
374 ndomains =
SIZE(subm)
375 DO idomain = 1, ndomains
376 CALL release_submatrix(subm(idomain))
379 CALL timestop(handle)
381 END SUBROUTINE release_submatrix_array
387 SUBROUTINE release_submatrix(subm)
389 TYPE(domain_submatrix_type),
INTENT(INOUT) :: subm
391 CHARACTER(len=*),
PARAMETER :: routinen =
'release_submatrix'
395 CALL timeset(routinen, handle)
405 IF (
ALLOCATED(subm%dbcsr_row))
THEN
406 DEALLOCATE (subm%dbcsr_row)
408 IF (
ALLOCATED(subm%dbcsr_col))
THEN
409 DEALLOCATE (subm%dbcsr_col)
411 IF (
ALLOCATED(subm%size_brow))
THEN
412 DEALLOCATE (subm%size_brow)
414 IF (
ALLOCATED(subm%size_bcol))
THEN
415 DEALLOCATE (subm%size_bcol)
417 IF (
ALLOCATED(subm%mdata))
THEN
418 DEALLOCATE (subm%mdata)
421 CALL timestop(handle)
423 END SUBROUTINE release_submatrix
436 SUBROUTINE multiply_submatrices_once(transA, transB, alpha, A, B, beta, C)
438 CHARACTER,
INTENT(IN) :: transa, transb
439 REAL(kind=
dp),
INTENT(IN) :: alpha
440 TYPE(domain_submatrix_type),
INTENT(IN) :: a, b
441 REAL(kind=
dp),
INTENT(IN) :: beta
442 TYPE(domain_submatrix_type),
INTENT(INOUT) :: c
444 CHARACTER(len=*),
PARAMETER :: routinen =
'multiply_submatrices_once'
446 INTEGER :: cs1, cs2, handle, icol, irow, k, k1, &
447 lda, ldb, ldc, m, mblocks, n, nblocks
448 LOGICAL :: nota, notb
450 CALL timeset(routinen, handle)
452 cpassert(a%domain .GT. 0)
453 cpassert(b%domain .GT. 0)
454 cpassert(c%domain .GT. 0)
456 lda =
SIZE(a%mdata, 1)
457 ldb =
SIZE(b%mdata, 1)
459 nota = (transa .EQ.
'N') .OR. (transa .EQ.
'n')
460 notb = (transb .EQ.
'N') .OR. (transb .EQ.
'n')
490 IF (
ALLOCATED(c%dbcsr_row))
THEN
491 DEALLOCATE (c%dbcsr_row)
493 ALLOCATE (c%dbcsr_row(c%nbrows))
494 IF (
ALLOCATED(c%dbcsr_col))
THEN
495 DEALLOCATE (c%dbcsr_col)
497 ALLOCATE (c%dbcsr_col(c%nbcols))
498 IF (
ALLOCATED(c%size_brow))
THEN
499 DEALLOCATE (c%size_brow)
501 ALLOCATE (c%size_brow(c%nbrows))
502 IF (
ALLOCATED(c%size_bcol))
THEN
503 DEALLOCATE (c%size_bcol)
505 ALLOCATE (c%size_bcol(c%nbcols))
507 DO irow = 1, c%nbrows
509 c%dbcsr_row(irow) = a%dbcsr_row(irow)
510 c%size_brow(irow) = a%size_brow(irow)
512 c%dbcsr_row(irow) = a%dbcsr_col(irow)
513 c%size_brow(irow) = a%size_bcol(irow)
517 DO icol = 1, c%nbcols
519 c%dbcsr_col(icol) = b%dbcsr_col(icol)
520 c%size_bcol(icol) = b%size_bcol(icol)
522 c%dbcsr_col(icol) = b%dbcsr_row(icol)
523 c%size_bcol(icol) = b%size_brow(icol)
527 IF (.NOT.
ALLOCATED(c%mdata))
THEN
529 cpassert(beta .EQ. 0.0_dp)
530 ALLOCATE (c%mdata(c%nrows, c%ncols))
532 cs1 =
SIZE(c%mdata, 1)
533 cs2 =
SIZE(c%mdata, 2)
534 IF ((c%nrows .NE. cs1) .OR. (c%ncols .NE. cs2))
THEN
536 cpassert(beta .EQ. 0.0_dp)
538 ALLOCATE (c%mdata(c%nrows, c%ncols))
544 CALL dgemm(transa, transb, m, n, k, alpha, a%mdata, lda, b%mdata, ldb, beta, c%mdata, ldc)
549 CALL timestop(handle)
551 END SUBROUTINE multiply_submatrices_once
563 SUBROUTINE multiply_submatrices_array(transA, transB, alpha, A, B, beta, C)
565 CHARACTER,
INTENT(IN) :: transa, transb
566 REAL(kind=
dp),
INTENT(IN) :: alpha
567 TYPE(domain_submatrix_type),
DIMENSION(:), &
569 REAL(kind=
dp),
INTENT(IN) :: beta
570 TYPE(domain_submatrix_type),
DIMENSION(:), &
573 CHARACTER(len=*),
PARAMETER :: routinen =
'multiply_submatrices_array'
575 INTEGER :: handle, idomain, idomaina, idomainb, &
576 ndomains, ndomainsb, ndomainsc
578 CALL timeset(routinen, handle)
584 cpassert(ndomains .EQ. ndomainsb)
585 cpassert(ndomainsb .EQ. ndomainsc)
587 DO idomain = 1, ndomains
589 idomaina = a(idomain)%domain
590 idomainb = b(idomain)%domain
592 cpassert(idomaina .EQ. idomainb)
594 c(idomain)%domain = idomaina
597 IF (idomaina .GT. 0)
THEN
598 CALL multiply_submatrices_once(transa, transb, alpha, a(idomain), b(idomain), beta, c(idomain))
603 CALL timestop(handle)
605 END SUBROUTINE multiply_submatrices_array
616 SUBROUTINE add_submatrices_once(alpha, A, beta, B, transB)
618 REAL(kind=
dp),
INTENT(IN) :: alpha
619 TYPE(domain_submatrix_type),
INTENT(INOUT) :: a
620 REAL(kind=
dp),
INTENT(IN) :: beta
621 TYPE(domain_submatrix_type),
INTENT(IN) :: b
622 CHARACTER,
INTENT(IN) :: transb
624 CHARACTER(len=*),
PARAMETER :: routinen =
'add_submatrices_once'
626 INTEGER :: c1, c2, handle, icol, r1, r2
629 CALL timeset(routinen, handle)
631 cpassert(a%domain .GT. 0)
632 cpassert(b%domain .GT. 0)
637 notb = (transb .EQ.
'N') .OR. (transb .EQ.
'n')
653 a%mdata(:, icol) = alpha*a%mdata(:, icol) + beta*b%mdata(:, icol)
657 a%mdata(:, icol) = alpha*a%mdata(:, icol) + beta*b%mdata(icol, :)
661 CALL timestop(handle)
663 END SUBROUTINE add_submatrices_once
673 SUBROUTINE add_submatrices_array(alpha, A, beta, B, transB)
675 REAL(kind=
dp),
INTENT(IN) :: alpha
676 TYPE(domain_submatrix_type),
DIMENSION(:), &
678 REAL(kind=
dp),
INTENT(IN) :: beta
679 TYPE(domain_submatrix_type),
DIMENSION(:), &
681 CHARACTER,
INTENT(IN) :: transb
683 CHARACTER(len=*),
PARAMETER :: routinen =
'add_submatrices_array'
685 INTEGER :: handle, idomain, idomaina, idomainb, &
688 CALL timeset(routinen, handle)
693 cpassert(ndomains .EQ. ndomainsb)
695 DO idomain = 1, ndomains
697 idomaina = a(idomain)%domain
698 idomainb = b(idomain)%domain
700 cpassert(idomaina .EQ. idomainb)
703 IF (idomaina .GT. 0)
THEN
704 CALL add_submatrices_once(alpha, a(idomain), beta, b(idomain), transb)
709 CALL timestop(handle)
711 END SUBROUTINE add_submatrices_array
723 TYPE(domain_submatrix_type),
DIMENSION(:), &
724 INTENT(IN) :: submatrices
725 REAL(kind=
dp),
INTENT(OUT) :: norm
727 CHARACTER(len=*),
PARAMETER :: routinen =
'maxnorm_submatrices'
729 INTEGER :: handle, idomain, ndomains
730 REAL(kind=
dp) :: curr_norm, send_norm
731 REAL(kind=
dp),
ALLOCATABLE,
DIMENSION(:) :: recv_norm
733 CALL timeset(routinen, handle)
737 ndomains =
SIZE(submatrices)
739 DO idomain = 1, ndomains
742 IF (submatrices(idomain)%domain .GT. 0)
THEN
743 curr_norm = maxval(abs(submatrices(idomain)%mdata))
744 IF (curr_norm .GT. send_norm) send_norm = curr_norm
750 ALLOCATE (recv_norm(submatrices(1)%nnodes))
751 CALL submatrices(1)%group%allgather(send_norm, recv_norm)
753 norm = maxval(recv_norm)
755 DEALLOCATE (recv_norm)
757 CALL timestop(handle)
770 SUBROUTINE trace_submatrices(A, B, trace)
772 TYPE(domain_submatrix_type),
DIMENSION(:), &
774 REAL(kind=
dp),
INTENT(OUT) :: trace
776 CHARACTER(len=*),
PARAMETER :: routinen =
'trace_submatrices'
778 INTEGER :: domaina, domainb, handle, idomain, &
780 REAL(kind=
dp) :: curr_trace, send_trace
781 REAL(kind=
dp),
ALLOCATABLE,
DIMENSION(:) :: recv_trace
783 CALL timeset(routinen, handle)
790 cpassert(ndomainsa .EQ. ndomainsb)
792 DO idomain = 1, ndomainsa
794 domaina = a(idomain)%domain
795 domainb = b(idomain)%domain
797 cpassert(domaina .EQ. domainb)
800 IF (domaina .GT. 0)
THEN
802 cpassert(a(idomain)%nrows .EQ. b(idomain)%nrows)
803 cpassert(a(idomain)%ncols .EQ. b(idomain)%ncols)
805 curr_trace = sum(a(idomain)%mdata(:, :)*b(idomain)%mdata(:, :))
806 send_trace = send_trace + curr_trace
813 ALLOCATE (recv_trace(a(1)%nnodes))
814 CALL a(1)%group%allgather(send_trace, recv_trace)
816 trace = sum(recv_trace)
818 DEALLOCATE (recv_trace)
820 CALL timestop(handle)
822 END SUBROUTINE trace_submatrices
838 node_of_domain, job_type)
840 TYPE(dbcsr_type),
INTENT(IN) :: matrix
841 TYPE(domain_submatrix_type),
DIMENSION(:), &
842 INTENT(INOUT) :: submatrix
843 TYPE(dbcsr_type),
INTENT(IN) :: distr_pattern
844 TYPE(domain_map_type),
INTENT(IN) :: domain_map
845 INTEGER,
DIMENSION(:),
INTENT(IN) :: node_of_domain
846 INTEGER,
INTENT(IN) :: job_type
848 CHARACTER(len=*),
PARAMETER :: routinen =
'construct_submatrices'
850 CHARACTER :: matrix_type
851 INTEGER :: block_node, block_offset, col, col_offset, col_size, dest_node, groupid, handle, &
852 iblock, icol, idomain, index_col, index_ec, index_er, index_row, index_sc, index_sr, &
853 inode, ldesc, mynode, nblkcols_tot, nblkrows_tot, ndomains, ndomains2, nnodes, &
854 recv_size2_total, recv_size_total, row, row_size, send_size2_total, send_size_total, &
855 smcol, smrow, start_data
856 INTEGER,
ALLOCATABLE,
DIMENSION(:) :: first_col, first_row, offset2_block, offset_block, &
857 recv_data2, recv_offset2_cpu, recv_offset_cpu, recv_size2_cpu, recv_size_cpu, send_data2, &
858 send_offset2_cpu, send_offset_cpu, send_size2_cpu, send_size_cpu
859 INTEGER,
ALLOCATABLE,
DIMENSION(:, :) :: recv_descriptor, send_descriptor
860 INTEGER,
DIMENSION(:),
POINTER :: col_blk_size, row_blk_size
861 LOGICAL :: found, transp
862 REAL(kind=
dp) :: antifactor
863 REAL(kind=
dp),
ALLOCATABLE,
DIMENSION(:) :: recv_data, send_data
864 REAL(kind=
dp),
DIMENSION(:),
POINTER :: block_p
865 TYPE(dbcsr_distribution_type) :: pattern_dist
866 TYPE(mp_comm_type) :: group
873 CALL timeset(routinen, handle)
875 nblkrows_tot = dbcsr_nblkrows_total(matrix)
876 nblkcols_tot = dbcsr_nblkcols_total(matrix)
878 ndomains = nblkcols_tot
879 CALL dbcsr_get_info(distr_pattern, distribution=pattern_dist)
880 CALL dbcsr_distribution_get(pattern_dist, numnodes=nnodes, group=groupid, mynode=mynode)
882 CALL group%set_handle(groupid)
884 matrix_type = dbcsr_get_matrix_type(matrix)
887 ALLOCATE (send_descriptor(ldesc, nnodes))
888 ALLOCATE (recv_descriptor(ldesc, nnodes))
889 send_descriptor(:, :) = 0
893 DO idomain = 1, ndomains
895 dest_node = node_of_domain(idomain)
899 IF (idomain .GT. 1) index_sr = domain_map%index1(idomain - 1)
900 index_er = domain_map%index1(idomain) - 1
902 DO index_row = index_sr, index_er
904 row = domain_map%pairs(index_row, 1)
909 IF (idomain .GT. 1) index_sc = domain_map%index1(idomain - 1)
910 index_ec = domain_map%index1(idomain) - 1
917 DO index_col = index_sc, index_ec
920 col = domain_map%pairs(index_col, 1)
926 CALL dbcsr_get_stored_coordinates(matrix, &
927 row, col, block_node)
928 IF (block_node .EQ. mynode)
THEN
929 CALL dbcsr_get_block_p(matrix, row, col, block_p, found, row_size, col_size)
931 send_descriptor(1, dest_node + 1) = send_descriptor(1, dest_node + 1) + 1
932 send_descriptor(2, dest_node + 1) = send_descriptor(2, dest_node + 1) + &
971 CALL group%alltoall(send_descriptor, recv_descriptor, ldesc)
973 ALLOCATE (send_size_cpu(nnodes), send_offset_cpu(nnodes))
974 send_offset_cpu(1) = 0
975 send_size_cpu(1) = send_descriptor(2, 1)
977 send_size_cpu(inode) = send_descriptor(2, inode)
978 send_offset_cpu(inode) = send_offset_cpu(inode - 1) + &
979 send_size_cpu(inode - 1)
981 send_size_total = send_offset_cpu(nnodes) + send_size_cpu(nnodes)
983 ALLOCATE (recv_size_cpu(nnodes), recv_offset_cpu(nnodes))
984 recv_offset_cpu(1) = 0
985 recv_size_cpu(1) = recv_descriptor(2, 1)
987 recv_size_cpu(inode) = recv_descriptor(2, inode)
988 recv_offset_cpu(inode) = recv_offset_cpu(inode - 1) + &
989 recv_size_cpu(inode - 1)
991 recv_size_total = recv_offset_cpu(nnodes) + recv_size_cpu(nnodes)
993 ALLOCATE (send_size2_cpu(nnodes), send_offset2_cpu(nnodes))
994 send_offset2_cpu(1) = 0
995 send_size2_cpu(1) = 2*send_descriptor(1, 1)
997 send_size2_cpu(inode) = 2*send_descriptor(1, inode)
998 send_offset2_cpu(inode) = send_offset2_cpu(inode - 1) + &
999 send_size2_cpu(inode - 1)
1001 send_size2_total = send_offset2_cpu(nnodes) + send_size2_cpu(nnodes)
1003 ALLOCATE (recv_size2_cpu(nnodes), recv_offset2_cpu(nnodes))
1004 recv_offset2_cpu(1) = 0
1005 recv_size2_cpu(1) = 2*recv_descriptor(1, 1)
1006 DO inode = 2, nnodes
1007 recv_size2_cpu(inode) = 2*recv_descriptor(1, inode)
1008 recv_offset2_cpu(inode) = recv_offset2_cpu(inode - 1) + &
1009 recv_size2_cpu(inode - 1)
1011 recv_size2_total = recv_offset2_cpu(nnodes) + recv_size2_cpu(nnodes)
1013 DEALLOCATE (send_descriptor)
1014 DEALLOCATE (recv_descriptor)
1017 ALLOCATE (send_data(send_size_total))
1018 ALLOCATE (recv_data(recv_size_total))
1019 ALLOCATE (send_data2(send_size2_total))
1020 ALLOCATE (recv_data2(recv_size2_total))
1021 ALLOCATE (offset_block(nnodes))
1022 ALLOCATE (offset2_block(nnodes))
1024 offset2_block(:) = 0
1026 DO idomain = 1, ndomains
1028 dest_node = node_of_domain(idomain)
1032 IF (idomain .GT. 1) index_sr = domain_map%index1(idomain - 1)
1033 index_er = domain_map%index1(idomain) - 1
1035 DO index_row = index_sr, index_er
1037 row = domain_map%pairs(index_row, 1)
1042 IF (idomain .GT. 1) index_sc = domain_map%index1(idomain - 1)
1043 index_ec = domain_map%index1(idomain) - 1
1050 DO index_col = index_sc, index_ec
1053 col = domain_map%pairs(index_col, 1)
1059 CALL dbcsr_get_stored_coordinates(matrix, &
1060 row, col, block_node)
1061 IF (block_node .EQ. mynode)
THEN
1062 CALL dbcsr_get_block_p(matrix, row, col, block_p, found, row_size, col_size)
1073 col_offset = row_size*col_size
1074 start_data = send_offset_cpu(dest_node + 1) + &
1075 offset_block(dest_node + 1)
1076 send_data(start_data + 1:start_data + col_offset) = &
1077 block_p(1:col_offset)
1078 offset_block(dest_node + 1) = offset_block(dest_node + 1) + col_offset
1080 send_data2(send_offset2_cpu(dest_node + 1) + &
1081 offset2_block(dest_node + 1) + 1) = row
1082 send_data2(send_offset2_cpu(dest_node + 1) + &
1083 offset2_block(dest_node + 1) + 2) = col
1084 offset2_block(dest_node + 1) = offset2_block(dest_node + 1) + 2
1134 CALL group%alltoall(send_data, send_size_cpu, send_offset_cpu, &
1135 recv_data, recv_size_cpu, recv_offset_cpu)
1137 CALL group%alltoall(send_data2, send_size2_cpu, send_offset2_cpu, &
1138 recv_data2, recv_size2_cpu, recv_offset2_cpu)
1140 DEALLOCATE (send_size_cpu, send_offset_cpu)
1141 DEALLOCATE (send_size2_cpu, send_offset2_cpu)
1142 DEALLOCATE (send_data)
1143 DEALLOCATE (send_data2)
1144 DEALLOCATE (offset_block)
1145 DEALLOCATE (offset2_block)
1148 CALL dbcsr_get_info(matrix, col_blk_size=col_blk_size, row_blk_size=row_blk_size)
1152 ndomains2 =
SIZE(submatrix)
1153 IF (ndomains2 .NE. ndomains)
THEN
1154 cpabort(
"wrong submatrix size")
1156 CALL release_submatrices(submatrix)
1157 submatrix(:)%nnodes = nnodes
1158 submatrix(:)%group = group
1159 submatrix(:)%nrows = 0
1160 submatrix(:)%ncols = 0
1162 ALLOCATE (first_row(nblkrows_tot), first_col(nblkcols_tot))
1163 submatrix(:)%domain = -1
1164 DO idomain = 1, ndomains
1165 dest_node = node_of_domain(idomain)
1169 IF (dest_node .EQ. mynode)
THEN
1170 submatrix(idomain)%domain = idomain
1171 submatrix(idomain)%nbrows = 0
1172 submatrix(idomain)%nbcols = 0
1177 IF (idomain .GT. 1) index_sr = domain_map%index1(idomain - 1)
1178 index_er = domain_map%index1(idomain) - 1
1179 DO index_row = index_sr, index_er
1180 row = domain_map%pairs(index_row, 1)
1183 first_row(row) = submatrix(idomain)%nrows + 1
1184 submatrix(idomain)%nrows = submatrix(idomain)%nrows + row_blk_size(row)
1185 submatrix(idomain)%nbrows = submatrix(idomain)%nbrows + 1
1188 ALLOCATE (submatrix(idomain)%dbcsr_row(submatrix(idomain)%nbrows))
1189 ALLOCATE (submatrix(idomain)%size_brow(submatrix(idomain)%nbrows))
1193 IF (idomain .GT. 1) index_sr = domain_map%index1(idomain - 1)
1194 index_er = domain_map%index1(idomain) - 1
1195 DO index_row = index_sr, index_er
1196 row = domain_map%pairs(index_row, 1)
1199 submatrix(idomain)%dbcsr_row(smrow) = row
1200 submatrix(idomain)%size_brow(smrow) = row_blk_size(row)
1210 IF (idomain .GT. 1) index_sc = domain_map%index1(idomain - 1)
1211 index_ec = domain_map%index1(idomain) - 1
1217 DO index_col = index_sc, index_ec
1219 col = domain_map%pairs(index_col, 1)
1230 first_col(col) = submatrix(idomain)%ncols + 1
1231 submatrix(idomain)%ncols = submatrix(idomain)%ncols + col_blk_size(col)
1232 submatrix(idomain)%nbcols = submatrix(idomain)%nbcols + 1
1236 ALLOCATE (submatrix(idomain)%dbcsr_col(submatrix(idomain)%nbcols))
1237 ALLOCATE (submatrix(idomain)%size_bcol(submatrix(idomain)%nbcols))
1244 IF (idomain .GT. 1) index_sc = domain_map%index1(idomain - 1)
1245 index_ec = domain_map%index1(idomain) - 1
1251 DO index_col = index_sc, index_ec
1253 col = domain_map%pairs(index_col, 1)
1259 submatrix(idomain)%dbcsr_col(smcol) = col
1260 submatrix(idomain)%size_bcol(smcol) = col_blk_size(col)
1265 ALLOCATE (submatrix(idomain)%mdata( &
1266 submatrix(idomain)%nrows, &
1267 submatrix(idomain)%ncols))
1268 submatrix(idomain)%mdata(:, :) = 0.0_dp
1269 DO inode = 1, nnodes
1271 DO iblock = 1, recv_size2_cpu(inode)/2
1273 row = recv_data2(recv_offset2_cpu(inode) + (iblock - 1)*2 + 1)
1274 col = recv_data2(recv_offset2_cpu(inode) + (iblock - 1)*2 + 2)
1276 IF ((first_col(col) .NE. -1) .AND. (first_row(row) .NE. -1))
THEN
1278 start_data = recv_offset_cpu(inode) + block_offset + 1
1279 DO icol = 0, col_blk_size(col) - 1
1280 submatrix(idomain)%mdata(first_row(row): &
1281 first_row(row) + row_blk_size(row) - 1, &
1282 first_col(col) + icol) = &
1283 recv_data(start_data:start_data + row_blk_size(row) - 1)
1284 start_data = start_data + row_blk_size(row)
1287 IF (matrix_type == dbcsr_type_symmetric .OR. &
1288 matrix_type == dbcsr_type_antisymmetric)
THEN
1291 IF (matrix_type == dbcsr_type_antisymmetric)
THEN
1292 antifactor = -1.0_dp
1294 start_data = recv_offset_cpu(inode) + block_offset + 1
1295 DO icol = 0, col_blk_size(col) - 1
1296 submatrix(idomain)%mdata(first_row(col) + icol, &
1298 first_col(row) + row_blk_size(row) - 1) = &
1299 antifactor*recv_data(start_data: &
1300 start_data + row_blk_size(row) - 1)
1301 start_data = start_data + row_blk_size(row)
1303 ELSE IF (matrix_type == dbcsr_type_no_symmetry)
THEN
1305 cpabort(
"matrix type is NYI")
1309 block_offset = block_offset + col_blk_size(col)*row_blk_size(row)
1315 DEALLOCATE (recv_size_cpu, recv_offset_cpu)
1316 DEALLOCATE (recv_size2_cpu, recv_offset2_cpu)
1317 DEALLOCATE (recv_data)
1318 DEALLOCATE (recv_data2)
1320 DEALLOCATE (first_row, first_col)
1322 CALL timestop(handle)
1337 TYPE(dbcsr_type),
INTENT(INOUT) :: matrix
1338 TYPE(domain_submatrix_type),
DIMENSION(:), &
1339 INTENT(IN) :: submatrix
1340 TYPE(dbcsr_type),
INTENT(IN) :: distr_pattern
1342 CHARACTER(len=*),
PARAMETER :: routinen =
'construct_dbcsr_from_submatrices'
1344 CHARACTER :: matrix_type
1345 INTEGER :: block_offset, col, col_offset, colsize, dest_node, groupid, handle, iblock, icol, &
1346 idomain, inode, irow_subm, ldesc, mynode, nblkcols_tot, nblkrows_tot, ndomains, &
1347 ndomains2, nnodes, recv_size2_total, recv_size_total, row, rowsize, send_size2_total, &
1348 send_size_total, smroff, start_data, unit_nr
1349 INTEGER,
ALLOCATABLE,
DIMENSION(:) :: offset2_block, offset_block, recv_data2, &
1350 recv_offset2_cpu, recv_offset_cpu, recv_size2_cpu, recv_size_cpu, send_data2, &
1351 send_offset2_cpu, send_offset_cpu, send_size2_cpu, send_size_cpu
1352 INTEGER,
ALLOCATABLE,
DIMENSION(:, :) :: recv_descriptor, send_descriptor
1353 INTEGER,
DIMENSION(:),
POINTER :: col_blk_size, row_blk_size
1355 REAL(kind=
dp),
ALLOCATABLE,
DIMENSION(:) :: recv_data, send_data
1356 REAL(kind=
dp),
DIMENSION(:, :),
POINTER :: block_p
1357 TYPE(cp_logger_type),
POINTER :: logger
1358 TYPE(dbcsr_distribution_type) :: pattern_dist
1359 TYPE(dbcsr_iterator_type) :: iter
1360 TYPE(mp_comm_type) :: group
1362 CALL timeset(routinen, handle)
1366 IF (logger%para_env%is_source())
THEN
1372 nblkrows_tot = dbcsr_nblkrows_total(matrix)
1373 nblkcols_tot = dbcsr_nblkcols_total(matrix)
1375 ndomains = nblkcols_tot
1376 ndomains2 =
SIZE(submatrix)
1378 IF (ndomains .NE. ndomains2)
THEN
1379 cpabort(
"domain mismatch")
1382 CALL dbcsr_get_info(distr_pattern, distribution=pattern_dist)
1383 CALL dbcsr_distribution_get(pattern_dist, numnodes=nnodes, group=groupid, mynode=mynode)
1385 CALL group%set_handle(groupid)
1387 matrix_type = dbcsr_get_matrix_type(matrix)
1388 IF (matrix_type .NE. dbcsr_type_no_symmetry)
THEN
1389 cpabort(
"only non-symmetric matrices so far")
1393 CALL dbcsr_iterator_start(iter, matrix)
1394 DO WHILE (dbcsr_iterator_blocks_left(iter))
1395 CALL dbcsr_iterator_next_block(iter, row, col, block_p)
1396 block_p(:, :) = 0.0_dp
1398 CALL dbcsr_iterator_stop(iter)
1399 CALL dbcsr_filter(matrix, 0.1_dp)
1401 CALL dbcsr_work_create(matrix, work_mutable=.true.)
1404 ALLOCATE (send_descriptor(ldesc, nnodes))
1405 ALLOCATE (recv_descriptor(ldesc, nnodes))
1406 send_descriptor(:, :) = 0
1409 DO idomain = 1, ndomains
1411 IF (submatrix(idomain)%domain .GT. 0)
THEN
1413 DO irow_subm = 1, submatrix(idomain)%nbrows
1415 IF (submatrix(idomain)%nbcols .NE. 1)
THEN
1416 cpabort(
"corrupt submatrix structure")
1419 row = submatrix(idomain)%dbcsr_row(irow_subm)
1420 col = submatrix(idomain)%dbcsr_col(1)
1422 IF (col .NE. idomain)
THEN
1423 cpabort(
"corrupt submatrix structure")
1427 CALL dbcsr_get_stored_coordinates(distr_pattern, &
1428 row, idomain, dest_node)
1430 send_descriptor(1, dest_node + 1) = send_descriptor(1, dest_node + 1) + 1
1431 send_descriptor(2, dest_node + 1) = send_descriptor(2, dest_node + 1) + &
1432 submatrix(idomain)%size_brow(irow_subm)* &
1433 submatrix(idomain)%size_bcol(1)
1442 CALL group%alltoall(send_descriptor, recv_descriptor, ldesc)
1444 ALLOCATE (send_size_cpu(nnodes), send_offset_cpu(nnodes))
1445 send_offset_cpu(1) = 0
1446 send_size_cpu(1) = send_descriptor(2, 1)
1447 DO inode = 2, nnodes
1448 send_size_cpu(inode) = send_descriptor(2, inode)
1449 send_offset_cpu(inode) = send_offset_cpu(inode - 1) + &
1450 send_size_cpu(inode - 1)
1452 send_size_total = send_offset_cpu(nnodes) + send_size_cpu(nnodes)
1454 ALLOCATE (recv_size_cpu(nnodes), recv_offset_cpu(nnodes))
1455 recv_offset_cpu(1) = 0
1456 recv_size_cpu(1) = recv_descriptor(2, 1)
1457 DO inode = 2, nnodes
1458 recv_size_cpu(inode) = recv_descriptor(2, inode)
1459 recv_offset_cpu(inode) = recv_offset_cpu(inode - 1) + &
1460 recv_size_cpu(inode - 1)
1462 recv_size_total = recv_offset_cpu(nnodes) + recv_size_cpu(nnodes)
1464 ALLOCATE (send_size2_cpu(nnodes), send_offset2_cpu(nnodes))
1465 send_offset2_cpu(1) = 0
1466 send_size2_cpu(1) = 2*send_descriptor(1, 1)
1467 DO inode = 2, nnodes
1468 send_size2_cpu(inode) = 2*send_descriptor(1, inode)
1469 send_offset2_cpu(inode) = send_offset2_cpu(inode - 1) + &
1470 send_size2_cpu(inode - 1)
1472 send_size2_total = send_offset2_cpu(nnodes) + send_size2_cpu(nnodes)
1474 ALLOCATE (recv_size2_cpu(nnodes), recv_offset2_cpu(nnodes))
1475 recv_offset2_cpu(1) = 0
1476 recv_size2_cpu(1) = 2*recv_descriptor(1, 1)
1477 DO inode = 2, nnodes
1478 recv_size2_cpu(inode) = 2*recv_descriptor(1, inode)
1479 recv_offset2_cpu(inode) = recv_offset2_cpu(inode - 1) + &
1480 recv_size2_cpu(inode - 1)
1482 recv_size2_total = recv_offset2_cpu(nnodes) + recv_size2_cpu(nnodes)
1484 DEALLOCATE (send_descriptor)
1485 DEALLOCATE (recv_descriptor)
1488 ALLOCATE (send_data(send_size_total))
1489 ALLOCATE (recv_data(recv_size_total))
1490 ALLOCATE (send_data2(send_size2_total))
1491 ALLOCATE (recv_data2(recv_size2_total))
1492 ALLOCATE (offset_block(nnodes))
1493 ALLOCATE (offset2_block(nnodes))
1495 offset2_block(:) = 0
1497 DO idomain = 1, ndomains
1499 IF (submatrix(idomain)%domain .GT. 0)
THEN
1502 DO irow_subm = 1, submatrix(idomain)%nbrows
1504 row = submatrix(idomain)%dbcsr_row(irow_subm)
1505 col = submatrix(idomain)%dbcsr_col(1)
1507 rowsize = submatrix(idomain)%size_brow(irow_subm)
1508 colsize = submatrix(idomain)%size_bcol(1)
1511 CALL dbcsr_get_stored_coordinates(distr_pattern, &
1512 row, idomain, dest_node)
1516 DO icol = 1, colsize
1517 start_data = send_offset_cpu(dest_node + 1) + &
1518 offset_block(dest_node + 1) + &
1520 send_data(start_data + 1:start_data + rowsize) = &
1521 submatrix(idomain)%mdata(smroff + 1:smroff + rowsize, icol)
1522 col_offset = col_offset + rowsize
1524 offset_block(dest_node + 1) = offset_block(dest_node + 1) + &
1527 send_data2(send_offset2_cpu(dest_node + 1) + &
1528 offset2_block(dest_node + 1) + 1) = row
1529 send_data2(send_offset2_cpu(dest_node + 1) + &
1530 offset2_block(dest_node + 1) + 2) = col
1531 offset2_block(dest_node + 1) = offset2_block(dest_node + 1) + 2
1533 smroff = smroff + rowsize
1542 CALL group%alltoall(send_data, send_size_cpu, send_offset_cpu, &
1543 recv_data, recv_size_cpu, recv_offset_cpu)
1545 CALL group%alltoall(send_data2, send_size2_cpu, send_offset2_cpu, &
1546 recv_data2, recv_size2_cpu, recv_offset2_cpu)
1548 DEALLOCATE (send_size_cpu, send_offset_cpu)
1549 DEALLOCATE (send_size2_cpu, send_offset2_cpu)
1550 DEALLOCATE (send_data)
1551 DEALLOCATE (send_data2)
1552 DEALLOCATE (offset_block)
1553 DEALLOCATE (offset2_block)
1556 CALL dbcsr_get_info(matrix, col_blk_size=col_blk_size, row_blk_size=row_blk_size)
1557 DO inode = 1, nnodes
1559 DO iblock = 1, recv_size2_cpu(inode)/2
1561 row = recv_data2(recv_offset2_cpu(inode) + (iblock - 1)*2 + 1)
1562 col = recv_data2(recv_offset2_cpu(inode) + (iblock - 1)*2 + 2)
1566 CALL dbcsr_reserve_block2d(matrix, row, col, block_p)
1567 cpassert(
ASSOCIATED(block_p))
1570 start_data = recv_offset_cpu(inode) + block_offset + 1
1571 DO icol = 1, col_blk_size(col)
1572 block_p(:, icol) = &
1573 recv_data(start_data:start_data + row_blk_size(row) - 1)
1574 start_data = start_data + row_blk_size(row)
1576 block_offset = block_offset + col_blk_size(col)*row_blk_size(row)
1580 DEALLOCATE (recv_size_cpu, recv_offset_cpu)
1581 DEALLOCATE (recv_size2_cpu, recv_offset2_cpu)
1582 DEALLOCATE (recv_data)
1583 DEALLOCATE (recv_data2)
1585 CALL dbcsr_finalize(matrix)
1587 CALL timestop(handle)
1598 TYPE(domain_submatrix_type),
DIMENSION(:), &
1599 INTENT(IN) :: submatrices
1600 TYPE(mp_comm_type),
INTENT(IN) :: mpgroup
1602 CHARACTER(len=*),
PARAMETER :: routinen =
'print_submatrices'
1604 CHARACTER(len=30) :: colstr, formatstr
1605 INTEGER :: handle, i, irow, n, ncols, nrows
1607 CALL timeset(routinen, handle)
1609 n =
SIZE(submatrices)
1612 nrows =
SIZE(submatrices(i)%mdata, 1)
1613 ncols =
SIZE(submatrices(i)%mdata, 2)
1614 WRITE (colstr, *) ncols
1615 formatstr =
'('//trim(adjustl(colstr))//
'F16.9)'
1616 IF (submatrices(i)%domain .GT. 0)
THEN
1617 WRITE (*, *)
"SUBMATRIX: ", i, nrows,
'x', ncols
1618 nrows =
SIZE(submatrices(i)%mdata, 1)
1620 WRITE (*, formatstr) submatrices(i)%mdata(irow, :)
1626 CALL timestop(handle)
1640 FUNCTION qblk_exists(map, row, col)
1642 TYPE(domain_map_type),
INTENT(IN) :: map
1643 INTEGER,
INTENT(IN) :: row, col
1644 LOGICAL :: qblk_exists
1646 INTEGER :: first, last, mid, ndomains
1650 ndomains =
SIZE(map%index1)
1652 qblk_exists = .false.
1653 IF (col .LT. 1 .OR. col .GT. ndomains)
RETURN
1655 IF (col .GT. 1) first = map%index1(col - 1)
1656 last = map%index1(col) - 1
1659 DO WHILE (last .GE. first)
1660 mid = first + (last - first)/2
1661 IF (map%pairs(mid, 1) .GT. row)
THEN
1663 ELSE IF (map%pairs(mid, 1) .LT. row)
THEN
1666 qblk_exists = .true.
1675 END FUNCTION qblk_exists
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.
various routines to log and control the output. The idea is that decisions about where to log should ...
recursive integer function, public cp_logger_get_default_unit_nr(logger, local, skip_not_ionode)
asks the default unit number of the given logger. try to use cp_logger_get_unit_nr
type(cp_logger_type) function, pointer, public cp_get_default_logger()
returns the default logger
Subroutines to handle submatrices.
subroutine, public copy_submatrix_data(array, copy)
...
subroutine, public print_submatrices(submatrices, mpgroup)
...
subroutine, public construct_dbcsr_from_submatrices(matrix, submatrix, distr_pattern)
Constructs a DBCSR matrix from submatrices.
subroutine, public maxnorm_submatrices(submatrices, norm)
Computes the max norm of the collection of submatrices.
subroutine, public construct_submatrices(matrix, submatrix, distr_pattern, domain_map, node_of_domain, job_type)
Constructs submatrices for each ALMO domain by collecting distributed DBCSR blocks to local arrays.
Types to handle submatrices.
integer, parameter, public select_row_col
Defines the basic variable types.
integer, parameter, public dp
Interface to the message passing library MPI.
type(mp_comm_type), parameter, public mp_comm_null