compute a norm of a dbcsr matrix
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
type(dbcsr_type), | intent(inout) | :: | matrix |
the matrix |
FUNCTION dbcsr_gershgorin_norm(matrix) RESULT(norm)
!! compute a norm of a dbcsr matrix
TYPE(dbcsr_type), INTENT(INOUT) :: matrix
!! the matrix
REAL(KIND=real_8) :: norm
CHARACTER(len=*), PARAMETER :: routineN = 'dbcsr_gershgorin_norm'
COMPLEX(KIND=real_4), DIMENSION(:, :), POINTER :: data_c
COMPLEX(KIND=real_8), DIMENSION(:, :), POINTER :: data_z
INTEGER :: blk, col, col_offset, handle, i, j, nc, &
nr, row, row_offset
LOGICAL :: any_sym, tr
REAL(KIND=real_4), DIMENSION(:, :), POINTER :: data_r
REAL(KIND=real_8), DIMENSION(:, :), POINTER :: data_d
REAL(real_8), ALLOCATABLE, DIMENSION(:) :: buff_d
TYPE(dbcsr_iterator) :: iter
CALL timeset(routineN, handle)
nr = dbcsr_nfullrows_total(matrix)
nc = dbcsr_nfullcols_total(matrix)
any_sym = dbcsr_get_matrix_type(matrix) .EQ. dbcsr_type_symmetric .OR. &
dbcsr_get_matrix_type(matrix) .EQ. dbcsr_type_antisymmetric
IF (nr .NE. nc) &
DBCSR_ABORT("not a square matrix")
norm = 0.0_dp
ALLOCATE (buff_d(nr))
buff_d = 0.0_dp
CALL dbcsr_iterator_start(iter, matrix)
DO WHILE (dbcsr_iterator_blocks_left(iter))
SELECT CASE (dbcsr_get_data_type(matrix))
CASE (dbcsr_type_real_4)
CALL dbcsr_iterator_next_block(iter, row, col, data_r, tr, blk, &
row_offset=row_offset, col_offset=col_offset)
DO j = 1, SIZE(data_r, 2)
DO i = 1, SIZE(data_r, 1)
buff_d(row_offset + i - 1) = buff_d(row_offset + i - 1) + ABS(data_r(i, j))
IF (any_sym .AND. row .NE. col) &
buff_d(col_offset + j - 1) = buff_d(col_offset + j - 1) + ABS(data_r(i, j))
END DO
END DO
CASE (dbcsr_type_real_8)
CALL dbcsr_iterator_next_block(iter, row, col, data_d, tr, blk, &
row_offset=row_offset, col_offset=col_offset)
DO j = 1, SIZE(data_d, 2)
DO i = 1, SIZE(data_d, 1)
buff_d(row_offset + i - 1) = buff_d(row_offset + i - 1) + ABS(data_d(i, j))
IF (any_sym .AND. row .NE. col) &
buff_d(col_offset + j - 1) = buff_d(col_offset + j - 1) + ABS(data_d(i, j))
END DO
END DO
CASE (dbcsr_type_complex_4)
CALL dbcsr_iterator_next_block(iter, row, col, data_c, tr, blk, &
row_offset=row_offset, col_offset=col_offset)
DO j = 1, SIZE(data_c, 2)
DO i = 1, SIZE(data_c, 1)
buff_d(row_offset + i - 1) = buff_d(row_offset + i - 1) + ABS(data_c(i, j))
IF (any_sym .AND. row .NE. col) &
DBCSR_ABORT("Only nonsymmetric matrix so far")
! buff_d(col_offset+j-1) = buff_d(col_offset+j-1) + ABS(data_c(i,j))
END DO
END DO
CASE (dbcsr_type_complex_8)
CALL dbcsr_iterator_next_block(iter, row, col, data_z, tr, blk, &
row_offset=row_offset, col_offset=col_offset)
DO j = 1, SIZE(data_z, 2)
DO i = 1, SIZE(data_z, 1)
buff_d(row_offset + i - 1) = buff_d(row_offset + i - 1) + ABS(data_z(i, j))
IF (any_sym .AND. row .NE. col) &
DBCSR_ABORT("Only nonsymmetric matrix so far")
! buff_d(col_offset+j-1) = buff_d(col_offset+j-1) + ABS(data_z(i,j))
END DO
END DO
CASE DEFAULT
DBCSR_ABORT("Wrong data type")
END SELECT
END DO
CALL dbcsr_iterator_stop(iter)
CALL mp_sum(buff_d, dbcsr_mp_group(dbcsr_distribution_mp(matrix%dist)))
norm = MAXVAL(buff_d)
DEALLOCATE (buff_d)
CALL timestop(handle)
END FUNCTION dbcsr_gershgorin_norm