dbcsr_trace_s Subroutine

private subroutine dbcsr_trace_s(matrix_a, trace)

traces a DBCSR matrix

Arguments

TypeIntentOptionalAttributesName
type(dbcsr_type), intent(in) :: matrix_a

DBCSR matrix

real(kind=real_4), intent(inout) :: trace

the trace of the matrix


Contents

Source Code


Source Code

      SUBROUTINE dbcsr_trace_s (matrix_a, trace)
      !! traces a DBCSR matrix

         TYPE(dbcsr_type), INTENT(IN)               :: matrix_a
         !! DBCSR matrix
         REAL(kind=real_4), INTENT(INOUT)                   :: trace
         !! the trace of the matrix

         CHARACTER(len=*), PARAMETER :: routineN = 'dbcsr_trace_s'

         INTEGER                                  :: a_blk, a_col, a_col_size, &
                                                     a_nze, a_row, a_row_size, i, &
                                                     mynode, error_handle
         INTEGER, DIMENSION(:), POINTER           :: col_blk_size, row_blk_size, &
                                                     row_dist, col_dist
         REAL(kind=real_4), DIMENSION(:), POINTER           :: a_data, data_p
         INTEGER, DIMENSION(:, :), POINTER         :: pgrid
         TYPE(dbcsr_distribution_obj)             :: dist

!   ---------------------------------------------------------------------------
         CALL timeset(routineN, error_handle)

         row_blk_size => array_data(matrix_a%row_blk_size)
         col_blk_size => array_data(matrix_a%col_blk_size)
         IF (dbcsr_get_data_type(matrix_a) /= dbcsr_type_real_4) &
            DBCSR_ABORT("Incompatible data types")
         CALL dbcsr_get_data(matrix_a%data_area, data_p)
         dist = dbcsr_distribution(matrix_a)
         mynode = dbcsr_mp_mynode(dbcsr_distribution_mp(dist))
         pgrid => dbcsr_mp_pgrid(dbcsr_distribution_mp(dist))
         row_dist => dbcsr_distribution_row_dist(dist)
         col_dist => dbcsr_distribution_col_dist(dist)
         !
         ! let's go
         trace = REAL(0.0, real_4)
         DO a_row = 1, matrix_a%nblkrows_total
            a_row_size = row_blk_size(a_row)
            DO a_blk = matrix_a%row_p(a_row) + 1, matrix_a%row_p(a_row + 1)
               IF (a_blk .EQ. 0) CYCLE
               a_col = matrix_a%col_i(a_blk)
               IF (a_col .ne. a_row) CYCLE
               ! We must skip non-local blocks in a replicated matrix.
               IF (matrix_a%replication_type .NE. dbcsr_repl_full) THEN
                  IF (mynode .NE. checker_square_proc(a_row, a_col, pgrid, row_dist, col_dist)) &
                     CYCLE
               END IF
               a_col_size = col_blk_size(a_col)
               IF (a_row_size .NE. a_col_size) &
                  DBCSR_ABORT("is that a square matrix?")
               a_nze = a_row_size**2
               a_data => pointer_view(data_p, ABS(matrix_a%blk_p(a_blk)), &
                                      ABS(matrix_a%blk_p(a_blk)) + a_nze - 1)
               !data_a => matrix_a%data(ABS(matrix_a%blk_p(a_blk)):ABS(matrix_a%blk_p(a_blk))+a_nze-1)
               !
               ! let's trace the block
               DO i = 1, a_row_size
                  trace = trace + a_data((i - 1)*a_row_size + i)
               END DO
            END DO ! a_col
         END DO ! a_row
         !
         ! summe
         CALL mp_sum(trace, dbcsr_mp_group(dbcsr_distribution_mp(matrix_a%dist)))

         CALL timestop(error_handle)
      END SUBROUTINE dbcsr_trace_s