dbcsr_get_diag_d Subroutine

      SUBROUTINE dbcsr_get_diag_d (matrix, diag)
         TYPE(dbcsr_type), INTENT(IN)               :: matrix
         REAL(kind=real_8), DIMENSION(:), INTENT(OUT)         :: diag

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

         INTEGER                                            :: icol, irow, row_offset, handle, i
         LOGICAL                                            :: tr
         TYPE(dbcsr_iterator)                               :: iter
         REAL(kind=real_8), DIMENSION(:, :), POINTER                   :: block

         CALL timeset(routineN, handle)

         IF (dbcsr_get_data_type(matrix) /= dbcsr_type_real_8) &
            DBCSR_ABORT("Incompatible data types")

         IF (dbcsr_nfullrows_total(matrix) /= SIZE(diag)) &
            DBCSR_ABORT("Diagonal has wrong size")

         IF (.NOT. array_equality(matrix%row_blk_offset, matrix%col_blk_offset)) &
            DBCSR_ABORT("matrix not quadratic")

         diag(:) = 0.0_real_8

         CALL dbcsr_iterator_start(iter, matrix)
         DO WHILE (dbcsr_iterator_blocks_left(iter))
            CALL dbcsr_iterator_next_block(iter, irow, icol, block, tr, row_offset=row_offset)
            IF (irow /= icol) CYCLE

            IF (sIZE(block, 1) /= sIZE(block, 2)) &
               DBCSR_ABORT("Diagonal block non-squared")

            DO i = 1, sIZE(block, 1)
               diag(row_offset + i - 1) = block(i, i)
            END DO
         END DO
         CALL dbcsr_iterator_stop(iter)

         CALL timestop(handle)
      END SUBROUTINE dbcsr_get_diag_d