convert_brd_to_csr Subroutine

   SUBROUTINE convert_brd_to_csr(brd_mat, csr_mat)
      !! Convert a block-row distributed DBCSR matrix to a CSR matrix
      !! The DBCSR matrix must have a block structure consistent with the CSR matrix.

      TYPE(dbcsr_type), INTENT(IN)                       :: brd_mat
         !! block-row distributed DBCSR matrix
      TYPE(csr_type), INTENT(INOUT)                      :: csr_mat
         !! CSR matrix

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

      INTEGER :: blk, blkcol, blkrow, col_blk_offset, col_blk_size, csr_ind, data_type, dbcsr_ind, &
                 el_sum, handle, ind_blk_data, k, local_row_ind, m, n, nblkrows_total, node_row_offset, &
                 prev_blkrow, prev_row_blk_size, row_blk_offset, row_blk_size
      INTEGER, ALLOCATABLE, DIMENSION(:)                 :: nfullcol_blkrow
      INTEGER, DIMENSION(:), POINTER                     :: colind, csr_index, dbcsr_index, nzerow, &
                                                            rowptr
      LOGICAL                                            :: new_ind, tr
      TYPE(dbcsr_data_obj)                               :: block
      TYPE(dbcsr_iterator)                               :: iter

      CALL timeset(routineN, handle)
      local_row_ind = 0
      dbcsr_ind = 0
      node_row_offset = 0
      NULLIFY (rowptr, colind, dbcsr_index, csr_index)

      dbcsr_index => csr_mat%dbcsr_mapping%csr_to_brd_ind
      csr_index => csr_mat%dbcsr_mapping%brd_to_csr_ind

      ! CSR indices are not recalculated if indices are already defined
      new_ind = .NOT. (csr_mat%has_indices)

      IF (.NOT. csr_mat%has_mapping) &
         DBCSR_ABORT("DBCSR mapping of CSR matrix must be defined")

      ! Calculate mapping between CSR matrix and DBCSR matrix if not yet defined
      !IF (.NOT. csr_mat%has_mapping ) THEN
      !  CALL csr_get_dbcsr_mapping (brd_mat, dbcsr_index, csr_index, nze)
      !ENDIF

      CALL dbcsr_get_info(brd_mat, nblkrows_total=nblkrows_total)
      ALLOCATE (nfullcol_blkrow(nblkrows_total))

      ! iteration over blocks without touching data,
      ! in order to get number of non-zero full columns in each block row
      CALL dbcsr_iterator_start(iter, brd_mat, read_only=.TRUE.)
      blkrow = 0
      nfullcol_blkrow = 0 ! number of non-zero full columns in each block row
      data_type = dbcsr_get_data_type(brd_mat)

      DO WHILE (dbcsr_iterator_blocks_left(iter))
         CALL dbcsr_iterator_next_block(iter, blkrow, blkcol, blk, col_size=col_blk_size, &
                                        row_offset=row_blk_offset)
         nfullcol_blkrow(blkrow) = nfullcol_blkrow(blkrow) + col_blk_size
         IF (blk .EQ. 1) THEN
            node_row_offset = row_blk_offset
         END IF
      END DO

      CALL dbcsr_iterator_stop(iter)

      ! Copy data from BRD matrix to CSR matrix and calculate CSR indices
      prev_blkrow = 0
      prev_row_blk_size = 0
      el_sum = 0 ! number of elements above current block row
      colind => csr_mat%colind_local
      rowptr => csr_mat%rowptr_local
      nzerow => csr_mat%nzerow_local
      CALL dbcsr_data_init(block)
      CALL dbcsr_data_new(block, data_type)

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

      IF (new_ind) rowptr(:) = 0 ! initialize to 0 in order to check which rows are 0 at a later time
      DO WHILE (dbcsr_iterator_blocks_left(iter))
         CALL dbcsr_iterator_next_block(iter, blkrow, blkcol, block, tr, &
                                        col_size=col_blk_size, row_size=row_blk_size, row_offset=row_blk_offset, &
                                        col_offset=col_blk_offset)

         IF (tr) &
            DBCSR_ABORT("DBCSR block data must not be transposed")

         IF (blkrow > prev_blkrow) THEN ! new block row
            local_row_ind = row_blk_offset - node_row_offset ! actually: local row index - 1
            IF (prev_blkrow .GT. 0) THEN
               el_sum = el_sum + nfullcol_blkrow(prev_blkrow)*prev_row_blk_size
            END IF
            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 = csr_index(dbcsr_ind) ! get CSR index for current element
               IF (csr_ind .GT. 0) THEN ! is non-zero element if csr_ind > 0
                  IF (new_ind) THEN
                     ! Calculate CSR column index
                     colind(csr_ind) = col_blk_offset + n - 1
                     ! Calculate CSR row pointer
                     ! (not accounting for zero elements inside non-zero blocks)
                     IF (rowptr(local_row_ind + m) .LE. 0) rowptr(local_row_ind + m) = &
                        rowptr(local_row_ind + m) + el_sum + 1 + nfullcol_blkrow(blkrow)*(m - 1)
                  END IF
                  ind_blk_data = (m + row_blk_size*(n - 1)) ! index of data inside DBCSR blocks
                  SELECT CASE (csr_mat%nzval_local%data_type)
                  CASE (dbcsr_type_real_4)
                     csr_mat%nzval_local%r_sp(csr_ind) = block%d%r_sp(ind_blk_data)
                  CASE (dbcsr_type_real_8)
                     csr_mat%nzval_local%r_dp(csr_ind) = block%d%r_dp(ind_blk_data)
                  CASE (dbcsr_type_complex_4)
                     csr_mat%nzval_local%c_sp(csr_ind) = block%d%c_sp(ind_blk_data)
                  CASE (dbcsr_type_complex_8)
                     csr_mat%nzval_local%c_dp(csr_ind) = block%d%c_dp(ind_blk_data)
                  END SELECT
               ELSE ! is zero element if ind = -1
                  ! CSR row pointer has to be corrected if element is zero
                  ! (subtract 1 from all subsequent row pointers)
                  IF (new_ind) rowptr(local_row_ind + m + 1:) = rowptr(local_row_ind + m + 1:) - 1
               END IF
            END DO
         END DO
         prev_blkrow = blkrow
         prev_row_blk_size = row_blk_size
      END DO

      IF (new_ind) THEN
         ! Repeat previous row pointer for row pointers to zero rows
         IF (csr_mat%nrows_local .GT. 0) rowptr(1) = 1
         DO k = 1, csr_mat%nrows_local
            IF (rowptr(k) .LE. 0) rowptr(k) = rowptr(k - 1)
         END DO

         rowptr(csr_mat%nrows_local + 1) = csr_mat%nze_local + 1
      END IF

      CALL csr_create_nzerow(csr_mat, nzerow)

      CALL dbcsr_iterator_stop(iter)
      CALL dbcsr_data_clear_pointer(block)
      CALL dbcsr_data_release(block)

      IF (new_ind) csr_mat%has_indices = .TRUE.

      CALL timestop(handle)
   END SUBROUTINE convert_brd_to_csr