dbcsr_copy_submatrix Subroutine

   SUBROUTINE dbcsr_copy_submatrix(matrix_b, matrix_a, name, &
                                   block_row_bounds, block_column_bounds, &
                                   shallow_data)
      !! Copy a submatrix.

      TYPE(dbcsr_type), INTENT(INOUT)                    :: matrix_b
         !! target DBCSR matrix
      TYPE(dbcsr_type), INTENT(IN)                       :: matrix_a
         !! source DBCSR matrix
      CHARACTER(LEN=*), INTENT(IN), OPTIONAL             :: name
         !! name of the new matrix
      INTEGER, DIMENSION(2), INTENT(IN), OPTIONAL        :: block_row_bounds, block_column_bounds
         !! rows to extract (array of size 2 holding the lower and upper inclusive bounds)
         !! columns to extract (array of size 2 holding the lower and upper inclusive bounds)
      LOGICAL, INTENT(IN), OPTIONAL                      :: shallow_data
         !! shallow data copy

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

      INTEGER                                            :: blk_p, col, handle, nblocks, new_blk, &
                                                            old_blk, row
      INTEGER, ALLOCATABLE, DIMENSION(:)                 :: blkp_list, col_list, row_list
      LOGICAL                                            :: shallow, tr
      TYPE(dbcsr_data_obj)                               :: data_block
      TYPE(dbcsr_iterator)                               :: iter

!   ---------------------------------------------------------------------------

      CALL timeset(routineN, handle)
      IF (PRESENT(shallow_data)) THEN
         shallow = shallow_data
      ELSE
         shallow = .FALSE.
      END IF
      ! Verify assumptions.
      IF (PRESENT(block_row_bounds)) THEN
         IF (SIZE(block_row_bounds) /= 2) &
            DBCSR_ABORT("Size of bounds specifier must be 2")
      END IF
      IF (PRESENT(block_column_bounds)) THEN
         IF (SIZE(block_column_bounds) /= 2) &
            DBCSR_ABORT("Size of bounds specifier must be 2")
      END IF
      ! Setup target matrix
      CALL dbcsr_create(matrix_b, name=name, template=matrix_a)
      CALL dbcsr_finalize(matrix_b)
      IF (.NOT. shallow) THEN
         ! Non-shallow copy uses the standard iterator on the source and
         ! block put on the target.
!
!$OMP PARALLEL DEFAULT (NONE) &
!$OMP          PRIVATE (data_block, iter, row, col, tr) &
!$OMP          SHARED (matrix_a, matrix_b,&
!$OMP                  block_row_bounds, block_column_bounds)
         CALL dbcsr_work_create(matrix_b, work_mutable=.FALSE.)
         CALL dbcsr_data_init(data_block)
         CALL dbcsr_data_new(data_block, dbcsr_get_data_type(matrix_a))
         CALL dbcsr_iterator_start(iter, matrix_a, dynamic=.TRUE., &
                                   dynamic_byrows=.TRUE.)
         DO WHILE (dbcsr_iterator_blocks_left(iter))
            CALL dbcsr_iterator_next_block(iter, row, col, data_block, tr)
            ! Only keep the block if they are within the specified bounds.
            IF (PRESENT(block_row_bounds)) THEN
               IF (row .LT. block_row_bounds(1)) CYCLE
               IF (row .GT. block_row_bounds(2)) CYCLE
            END IF
            IF (PRESENT(block_column_bounds)) THEN
               IF (col .LT. block_column_bounds(1)) CYCLE
               IF (col .GT. block_column_bounds(2)) CYCLE
            END IF
            CALL dbcsr_put_block(matrix_b, row, col, data_block, transposed=tr)
         END DO
         CALL dbcsr_iterator_stop(iter)
         CALL dbcsr_data_clear_pointer(data_block)
         CALL dbcsr_data_release(data_block)
         CALL dbcsr_finalize(matrix_b)
!$OMP END PARALLEL
      ELSE
         ! For the shallow copy the source matrix data is referenced.
         CALL dbcsr_switch_data_area(matrix_b, matrix_a%data_area)
         nblocks = dbcsr_get_num_blocks(matrix_a) ! High estimate.
         ! Shallow copy goes through source's data blocks and inserts
         ! the only the ones corresponding to the submatrix specifier
         ! into the target. Block pointers must remain the same as in
         ! the source.
         ALLOCATE (row_list(nblocks), col_list(nblocks), blkp_list(nblocks))
         !
         CALL dbcsr_iterator_start(iter, matrix_a)
         new_blk = 1
         DO WHILE (dbcsr_iterator_blocks_left(iter))
            CALL dbcsr_iterator_next_block(iter, row, col, &
                                           blk=old_blk, blk_p=blk_p)
            ! Only keep the block if they are within the specified bounds.
            IF (PRESENT(block_row_bounds)) THEN
               IF (row .LT. block_row_bounds(1)) CYCLE
               IF (row .GT. block_row_bounds(2)) CYCLE
            END IF
            IF (PRESENT(block_column_bounds)) THEN
               IF (col .LT. block_column_bounds(1)) CYCLE
               IF (col .GT. block_column_bounds(2)) CYCLE
            END IF
            row_list(new_blk) = row
            col_list(new_blk) = col
            blkp_list(new_blk) = blk_p
            new_blk = new_blk + 1
         END DO
         new_blk = new_blk - 1
         CALL dbcsr_iterator_stop(iter)
         CALL dbcsr_reserve_blocks(matrix_b, row_list(1:new_blk), &
                                   col_list(1:new_blk), blkp_list(1:new_blk))
      END IF
      !
      CALL timestop(handle)
   END SUBROUTINE dbcsr_copy_submatrix