Duplicates data in symmetric matrix to make it normal (w.r.t. data structure. Can handle symmetric/antisymmetric/hermitian/skewhermitian matrices
| Type | Intent | Optional | Attributes | Name | ||
|---|---|---|---|---|---|---|
| type(dbcsr_type), | intent(in) | :: | sm |
input symmetric matrix |
||
| type(dbcsr_type), | intent(inout) | :: | desm |
desymmetrized matrix |
||
| logical, | intent(in), | optional | :: | untransposed_data |
SUBROUTINE dbcsr_desymmetrize_deep(sm, desm, untransposed_data)
!! Duplicates data in symmetric matrix to make it normal (w.r.t. data
!! structure. Can handle symmetric/antisymmetric/hermitian/skewhermitian
!! matrices
TYPE(dbcsr_type), INTENT(IN) :: sm
!! input symmetric matrix
TYPE(dbcsr_type), INTENT(INOUT) :: desm
!! desymmetrized matrix
LOGICAL, INTENT(IN), OPTIONAL :: untransposed_data
CHARACTER(len=*), PARAMETER :: routineN = 'dbcsr_desymmetrize_deep'
INTEGER, PARAMETER :: idata = 2, imeta = 1, metalen = 3
INTEGER :: blk, blk_l, blk_p, blk_ps, blks, col, col_size, dst_p, handle, mp_group, &
nsymmetries, numproc, nze, pcol, prow, row, row_size, src_p, stored_col, stored_row, &
symmetry_i
INTEGER, ALLOCATABLE, DIMENSION(:) :: rd_disp, recv_meta, rm_disp, sd_disp, &
sdp, send_meta, sm_disp, smp
INTEGER, ALLOCATABLE, DIMENSION(:, :) :: recv_count, send_count, &
total_recv_count, total_send_count
INTEGER, DIMENSION(:), POINTER :: col_blk_size, col_dist, row_blk_size, &
row_dist
INTEGER, DIMENSION(:, :), POINTER :: blacs2mpi
LOGICAL :: make_untr, tr
TYPE(dbcsr_data_obj) :: recv_data, send_data
TYPE(dbcsr_distribution_obj) :: target_dist
TYPE(dbcsr_iterator) :: iter
TYPE(dbcsr_mp_obj) :: mp_obj
! ---------------------------------------------------------------------------
CALL timeset(routineN, handle)
IF (.NOT. dbcsr_valid_index(sm)) &
DBCSR_ABORT("Input matrix is invalid.")
IF (dbcsr_valid_index(desm)) THEN
CALL dbcsr_release(desm)
END IF
CALL dbcsr_create(desm, sm, matrix_type="N")
IF (PRESENT(untransposed_data)) THEN
make_untr = untransposed_data
ELSE
make_untr = .FALSE.
END IF
nsymmetries = 1
IF (sm%symmetry) THEN
nsymmetries = 2
END IF
row_blk_size => array_data(sm%row_blk_size)
col_blk_size => array_data(sm%col_blk_size)
target_dist = sm%dist
row_dist => dbcsr_distribution_row_dist(target_dist)
col_dist => dbcsr_distribution_col_dist(target_dist)
mp_obj = dbcsr_distribution_mp(target_dist)
blacs2mpi => dbcsr_mp_pgrid(mp_obj)
numproc = dbcsr_mp_numnodes(mp_obj)
mp_group = dbcsr_mp_group(mp_obj)
!
ALLOCATE (send_count(2, 0:numproc - 1))
ALLOCATE (recv_count(2, 0:numproc - 1))
ALLOCATE (total_send_count(2, 0:numproc - 1))
ALLOCATE (total_recv_count(2, 0:numproc - 1))
ALLOCATE (sdp(0:numproc - 1))
ALLOCATE (sd_disp(0:numproc - 1))
ALLOCATE (smp(0:numproc - 1))
ALLOCATE (sm_disp(0:numproc - 1))
ALLOCATE (rd_disp(0:numproc - 1))
ALLOCATE (rm_disp(0:numproc - 1))
!
desm%negate_real = sm%negate_real
desm%negate_imaginary = sm%negate_imaginary
! Count initial sizes for sending.
send_count(:, :) = 0
CALL dbcsr_iterator_start(iter, sm)
DO WHILE (dbcsr_iterator_blocks_left(iter))
CALL dbcsr_iterator_next_block(iter, row, col, blk, &
row_size=row_size, col_size=col_size)
DO symmetry_i = 1, nsymmetries
IF (symmetry_i .EQ. 1) THEN
stored_row = row; stored_col = col
ELSE
IF (row .EQ. col) CYCLE
stored_row = col; stored_col = row
END IF
! Where do we send this block?
prow = row_dist(stored_row)
pcol = col_dist(stored_col)
dst_p = blacs2mpi(prow, pcol)
nze = row_size*col_size
send_count(imeta, dst_p) = send_count(imeta, dst_p) + 1
send_count(idata, dst_p) = send_count(idata, dst_p) + nze
END DO ! symmetry_i
END DO ! iter
CALL dbcsr_iterator_stop(iter)
!
CALL mp_alltoall(send_count, recv_count, 2, mp_group)
!
! Allocate data structures needed for data exchange.
CALL dbcsr_data_init(recv_data)
CALL dbcsr_data_new(recv_data, dbcsr_get_data_type(sm), &
SUM(recv_count(idata, :)))
ALLOCATE (recv_meta(metalen*SUM(recv_count(imeta, :))))
CALL dbcsr_data_init(send_data)
CALL dbcsr_data_new(send_data, dbcsr_get_data_type(sm), &
SUM(send_count(idata, :)))
ALLOCATE (send_meta(metalen*SUM(send_count(imeta, :))))
!
! Fill in the meta data structures and copy the data.
DO dst_p = 0, numproc - 1
total_send_count(imeta, dst_p) = send_count(imeta, dst_p)
total_send_count(idata, dst_p) = send_count(idata, dst_p)
total_recv_count(imeta, dst_p) = recv_count(imeta, dst_p)
total_recv_count(idata, dst_p) = recv_count(idata, dst_p)
END DO
sd_disp = -1; sm_disp = -1
rd_disp = -1; rm_disp = -1
sd_disp(0) = 1; sm_disp(0) = 1
rd_disp(0) = 1; rm_disp(0) = 1
DO dst_p = 1, numproc - 1
sm_disp(dst_p) = sm_disp(dst_p - 1) &
+ metalen*total_send_count(imeta, dst_p - 1)
sd_disp(dst_p) = sd_disp(dst_p - 1) &
+ total_send_count(idata, dst_p - 1)
rm_disp(dst_p) = rm_disp(dst_p - 1) &
+ metalen*total_recv_count(imeta, dst_p - 1)
rd_disp(dst_p) = rd_disp(dst_p - 1) &
+ total_recv_count(idata, dst_p - 1)
END DO
sdp(:) = sd_disp
smp(:) = sm_disp
!
CALL dbcsr_iterator_start(iter, sm)
DO WHILE (dbcsr_iterator_blocks_left(iter))
CALL dbcsr_iterator_next_block(iter, row, col, blk, blk_p=blk_p, &
row_size=row_size, col_size=col_size)
DO symmetry_i = 1, nsymmetries
IF (symmetry_i .EQ. 1) THEN
stored_row = row; stored_col = col; tr = .FALSE.
ELSE
IF (row .EQ. col) CYCLE
stored_row = col; stored_col = row; tr = .TRUE.
END IF
! Where do we send this block?
prow = row_dist(stored_row)
pcol = col_dist(stored_col)
dst_p = blacs2mpi(prow, pcol)
nze = row_size*col_size
send_meta(smp(dst_p)) = stored_row
send_meta(smp(dst_p) + 1) = stored_col
send_meta(smp(dst_p) + 2) = sdp(dst_p) - sd_disp(dst_p) + 1
IF (make_untr) THEN
CALL dbcsr_block_partial_copy(dst=send_data, &
dst_rs=row_size, dst_cs=col_size, &
dst_tr=symmetry_i .GT. 1, &
dst_r_lb=1, dst_c_lb=1, dst_offset=sdp(dst_p) - 1, &
nrow=row_size, ncol=col_size, &
src=sm%data_area, &
src_rs=row_size, src_cs=col_size, &
src_tr=blk_p .LT. 0, &
src_r_lb=1, src_c_lb=1, &
src_offset=ABS(blk_p) - 1)
IF (tr) THEN
IF (sm%negate_real) THEN
CALL dbcsr_block_real_neg(dst=send_data, row_size=row_size, &
col_size=col_size, lb=sdp(dst_p))
END IF
IF (sm%negate_imaginary) &
CALL dbcsr_block_conjg(dst=send_data, row_size=row_size, &
col_size=col_size, lb=sdp(dst_p))
END IF
ELSE
CALL dbcsr_data_copy(send_data, (/sdp(dst_p)/), (/nze/), &
sm%data_area, (/ABS(blk_p)/), (/nze/))
IF (tr) &
send_meta(smp(dst_p) + 2) = -send_meta(smp(dst_p) + 2)
END IF
smp(dst_p) = smp(dst_p) + metalen
sdp(dst_p) = sdp(dst_p) + nze
END DO ! symmetry_i
END DO ! iter
CALL dbcsr_iterator_stop(iter)
! Exchange the data and metadata structures.
CALL hybrid_alltoall_any(send_data, &
total_send_count(idata, :), sd_disp(:) - 1, &
recv_data, &
total_recv_count(idata, :), rd_disp(:) - 1, mp_obj)
CALL mp_alltoall(send_meta(:), metalen*total_send_count(imeta, :), sm_disp(:) - 1, &
recv_meta(:), metalen*total_recv_count(imeta, :), rm_disp(:) - 1, &
mp_group)
! Now fill in the data.
CALL dbcsr_work_create(desm, &
SUM(recv_count(imeta, :)), &
SUM(recv_count(idata, :)), n=1, &
work_mutable=.FALSE.)
! Switch data data area of the work matrix with the received data
! (avoids copying).
CALL dbcsr_data_hold(recv_data)
CALL dbcsr_data_release(desm%wms(1)%data_area)
desm%wms(1)%data_area = recv_data
!
blk_ps = 1
blks = 1
! WRITE(*,*)rm_disp
! WRITE(*,*)recv_count
DO src_p = 0, numproc - 1
IF (careful_mod) THEN
IF ((blks - 1)*3 /= rm_disp(src_p) - 1) &
DBCSR_ABORT("Count mismatch")
END IF
blks = (rm_disp(src_p) - 1)/metalen + 1
DO blk_l = 1, recv_count(imeta, src_p)
stored_row = recv_meta(rm_disp(src_p) + metalen*(blk_l - 1))
stored_col = recv_meta(rm_disp(src_p) + metalen*(blk_l - 1) + 1)
blk_p = recv_meta(rm_disp(src_p) + metalen*(blk_l - 1) + 2)
desm%wms(1)%row_i(blks) = stored_row
desm%wms(1)%col_i(blks) = stored_col
desm%wms(1)%blk_p(blks) = SIGN(ABS(blk_p) + (rd_disp(src_p) - 1), blk_p)
nze = row_blk_size(ABS(stored_row)) &
*col_blk_size(stored_col)
blk_ps = blk_ps + nze
blks = blks + 1
END DO
END DO
!
desm%wms(1)%lastblk = blks - 1
desm%wms(1)%datasize = blk_ps - 1
CALL dbcsr_finalize(desm)
DEALLOCATE (send_count)
DEALLOCATE (recv_count)
DEALLOCATE (sdp); DEALLOCATE (sd_disp)
DEALLOCATE (smp); DEALLOCATE (sm_disp)
DEALLOCATE (rd_disp)
DEALLOCATE (rm_disp)
CALL dbcsr_data_release(recv_data)
DEALLOCATE (recv_meta)
CALL dbcsr_data_release(send_data)
DEALLOCATE (send_meta)
CALL timestop(handle)
END SUBROUTINE dbcsr_desymmetrize_deep