csr_get_dbcsr_mapping Subroutine

   SUBROUTINE csr_get_dbcsr_mapping(brd_mat, dbcsr_index, csr_index, csr_nze_local, &
                                    csr_sparsity_brd)
      !! create the mapping information between a block-row distributed DBCSR
      !! matrix and the corresponding CSR matrix

      TYPE(dbcsr_type), INTENT(IN)                       :: brd_mat
         !! the block-row distributed DBCSR matrix
      INTEGER, DIMENSION(:), INTENT(OUT), POINTER        :: dbcsr_index, csr_index
         !! csr to dbcsr index mapping
         !! dbcsr to csr index mapping
      INTEGER, INTENT(OUT)                               :: csr_nze_local
         !! number of local non-zero elements
      TYPE(dbcsr_type), INTENT(IN)                       :: csr_sparsity_brd
         !! sparsity of CSR matrix represented in BRD format

      CHARACTER(LEN=*), PARAMETER :: routineN = 'csr_get_dbcsr_mapping'

      INTEGER :: blk, blkcol, blkrow, col_blk_size, csr_ind, data_type, dbcsr_ind, el_sum, &
                 fullcol_sum_blkrow, handle, l, m, n, nblkrows_total, nze, prev_blk, prev_blkrow, &
                 prev_row_blk_size, row_blk_offset, row_blk_size
      INTEGER, ALLOCATABLE, DIMENSION(:)                 :: csr_nze, nfullcol_blkrow
      INTEGER, DIMENSION(:), POINTER                     :: dbcsr_index_nozeroes
      LOGICAL                                            :: tr
      TYPE(dbcsr_iterator)                               :: iter

      CALL timeset(routineN, handle)

      m = 0
      dbcsr_ind = 0
      fullcol_sum_blkrow = 0
      NULLIFY (dbcsr_index, csr_index)

      CALL dbcsr_get_info(brd_mat, nblkrows_total=nblkrows_total)
      nze = dbcsr_get_nze(brd_mat)

      ALLOCATE (nfullcol_blkrow(nblkrows_total))
      ALLOCATE (dbcsr_index(nze))
      ALLOCATE (csr_index(nze))

      CALL dbcsr_iterator_start(iter, brd_mat, read_only=.TRUE.)
      nfullcol_blkrow = 0 ! number of non-zero full columns in each block row
      prev_blk = 0

      DO WHILE (dbcsr_iterator_blocks_left(iter))
         CALL dbcsr_iterator_next_block(iter, blkrow, blkcol, blk, transposed=tr, &
                                        col_size=col_blk_size)

         IF (blk /= prev_blk + 1) &
            DBCSR_ABORT("iterator is required to traverse the blocks in a row-wise fashion")

         prev_blk = blk

         nfullcol_blkrow(blkrow) = nfullcol_blkrow(blkrow) + col_blk_size
         IF (tr) &
            DBCSR_ABORT("DBCSR block data must not be transposed")
      END DO
      CALL dbcsr_iterator_stop(iter)

      el_sum = 0 ! number of elements above current block row

      prev_blkrow = 0 ! store number and size of previous block row
      prev_row_blk_size = 0

      CALL dbcsr_iterator_start(iter, brd_mat, read_only=.TRUE.)

      DO WHILE (dbcsr_iterator_blocks_left(iter))

         CALL dbcsr_iterator_next_block(iter, blkrow, blkcol, blk, transposed=tr, &
                                        row_size=row_blk_size, col_size=col_blk_size, row_offset=row_blk_offset)

         IF (blkrow .GT. prev_blkrow) THEN ! new block row
            IF (prev_blkrow .GT. 0) THEN
               el_sum = el_sum + nfullcol_blkrow(prev_blkrow)*prev_row_blk_size
            END IF

            ! number of non-zero full columns on the left of current block:
            fullcol_sum_blkrow = 0

            dbcsr_ind = el_sum
         END IF
         DO n = 1, col_blk_size !nr of columns
            DO m = 1, row_blk_size !nr of rows
               dbcsr_ind = dbcsr_ind + 1
               csr_ind = (m - 1)*nfullcol_blkrow(blkrow) + fullcol_sum_blkrow + n + el_sum
               dbcsr_index(csr_ind) = dbcsr_ind
               csr_index(dbcsr_ind) = csr_ind
            END DO
         END DO
         fullcol_sum_blkrow = fullcol_sum_blkrow + col_blk_size
         prev_blkrow = blkrow
         prev_row_blk_size = row_blk_size
      END DO
      CALL dbcsr_iterator_stop(iter)

      ! remove BRD zero elements from CSR format
      data_type = dbcsr_get_data_type(csr_sparsity_brd)
      ALLOCATE (csr_nze(nze))

      SELECT CASE (data_type)
      CASE (dbcsr_type_real_4)
         csr_nze(:) = INT(csr_sparsity_brd%data_area%d%r_sp(1:nze))
      CASE (dbcsr_type_real_8)
         csr_nze(:) = INT(csr_sparsity_brd%data_area%d%r_dp(1:nze))
      CASE DEFAULT
         DBCSR_ABORT("CSR sparsity matrix must have a real datatype")
      END SELECT

      IF (ANY(csr_nze .EQ. 0)) THEN
         ALLOCATE (dbcsr_index_nozeroes(SUM(csr_nze)))
         m = 0 ! csr index if zeroes are excluded from CSR data
         DO l = 1, nze ! csr index if zeroes are included in CSR data
            IF (csr_nze(dbcsr_index(l)) .EQ. 0) THEN
               csr_index(dbcsr_index(l)) = -1
            ELSE
               m = m + 1
               dbcsr_index_nozeroes(m) = dbcsr_index(l)
               csr_index(dbcsr_index(l)) = m
            END IF
         END DO
         DEALLOCATE (dbcsr_index)
         dbcsr_index => dbcsr_index_nozeroes
      END IF

      IF (ANY(csr_nze .EQ. 0)) THEN
         csr_nze_local = m
      ELSE
         csr_nze_local = nze
      END IF

      CALL timestop(handle)
   END SUBROUTINE csr_get_dbcsr_mapping