dbcsr_make_index_list Subroutine

public subroutine dbcsr_make_index_list(matrix, thread_redist)

Converts BCSR index into list indices (similar to work matrices)

Arguments

TypeIntentOptionalAttributesName
type(dbcsr_type), intent(inout) :: matrix

matrix for which to make canonical index

logical, intent(in) :: thread_redist

make a thread subdistribution


Contents

Source Code


Source Code

   SUBROUTINE dbcsr_make_index_list(matrix, thread_redist)
      !! Converts BCSR index into list indices (similar to work matrices)

      TYPE(dbcsr_type), INTENT(INOUT)                    :: matrix
         !! matrix for which to make canonical index
      LOGICAL, INTENT(IN)                                :: thread_redist
         !! make a thread subdistribution

      CHARACTER(len=*), PARAMETER :: routineN = 'dbcsr_make_index_list'
      LOGICAL, PARAMETER                                 :: dbg = .FALSE.

      INTEGER                                            :: blk, error_handle, ithread, my_cnt, &
                                                            nblks, nrows, nthreads
      INTEGER, ALLOCATABLE, DIMENSION(:, :)              :: blki
      INTEGER, DIMENSION(0)                              :: zero_len_array
      INTEGER, DIMENSION(:), POINTER                     :: global_cols, local_rows, td, thr_c

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

      CALL timeset(routineN, error_handle)
      IF (.NOT. ASSOCIATED(matrix%row_p)) &
         DBCSR_ABORT("The row index must be initialized.")
      IF (matrix%list_indexing) &
         DBCSR_ABORT("List index already exists?")
      IF (matrix%bcsc) &
         DBCSR_ABORT("Not support for BCSC yet.")
      IF (matrix%nblks /= SIZE(matrix%col_i)) &
         DBCSR_ABORT("Block count mismatch.")
      IF (matrix%nblks /= SIZE(matrix%blk_p)) &
         DBCSR_ABORT("Block count mismatch")
      !
      IF (matrix%local_indexing) THEN
         IF (SIZE(matrix%row_p) - 1 /= matrix%nblkrows_local) &
            DBCSR_ABORT("Local row index incorrectly sized.")
      ELSE
         IF (SIZE(matrix%row_p) - 1 /= matrix%nblkrows_total) &
            DBCSR_ABORT("Global row index incorrectly sized")
      END IF
      !
      matrix%list_indexing = .TRUE.
      !
      IF (matrix%local_indexing) THEN
         nrows = matrix%nblkrows_local
      ELSE
         nrows = matrix%nblkrows_total
      END IF
      !
      nblks = matrix%nblks
      ALLOCATE (blki(3, nblks))
      CALL dbcsr_expand_row_index_2d(matrix%row_p, blki, nrows, 1)
      IF (matrix%local_indexing) THEN
         global_cols => array_data(matrix%global_cols)
         ! If local indexing is enabled, then the rows but not the
         ! columns are already localized
         IF (dbg) THEN
            WRITE (*, *) routineN//" Making local columns"
            WRITE (*, '(10(1X,i7))') global_cols
            WRITE (*, *) 'local'
            WRITE (*, '(10(1X,i7))') array_data(matrix%local_cols)
         END IF
         DO blk = 1, nblks
            blki(2, blk) = global_cols(matrix%col_i(blk))
            blki(3, blk) = matrix%blk_p(blk)
         END DO
      ELSE
         IF (dbg) WRITE (*, *) routineN//" Not making local columns"
         DO blk = 1, nblks
            blki(2, blk) = matrix%col_i(blk)
            blki(3, blk) = matrix%blk_p(blk)
         END DO
      END IF
      CALL dbcsr_clearfrom_index_array(matrix, dbcsr_slot_row_p)
      CALL dbcsr_clearfrom_index_array(matrix, dbcsr_slot_col_i)
      CALL dbcsr_clearfrom_index_array(matrix, dbcsr_slot_blk_p)
      nthreads = 0
!$    nthreads = OMP_GET_MAX_THREADS()
      !
      ! Reshuffle according to threads
      IF (nthreads .GT. 0 .AND. thread_redist) THEN
         td => array_data(dbcsr_distribution_thread_dist(dbcsr_distribution(matrix)))
         IF (matrix%local_indexing) THEN
            local_rows => array_data(matrix%local_rows)
         END IF
!$OMP        PARALLEL DEFAULT (NONE) &
!$OMP        PRIVATE (my_cnt, ithread, blk) &
!$OMP        SHARED (td, blki, nthreads, thr_c, nblks, matrix,local_rows)
         !
         ithread = 0
!$       ithread = OMP_GET_THREAD_NUM()
!$OMP        MASTER
!$       nthreads = OMP_GET_NUM_THREADS()
         CALL dbcsr_addto_index_array(matrix, dbcsr_slot_thr_c, &
                                      reservation=nthreads + 1, extra=3*nblks)
         thr_c => matrix%thr_c
         CALL dbcsr_addto_index_array(matrix, dbcsr_slot_coo_l, &
                                      reservation=3*nblks)
!$OMP        END MASTER
!$OMP        BARRIER
         my_cnt = 0
         IF (matrix%local_indexing) THEN
            my_cnt = COUNT(td(local_rows(blki(1, :))) .EQ. ithread)
         ELSE
            my_cnt = COUNT(td(blki(1, :)) .EQ. ithread)
         END IF
         !DO blk = 1, nblks
         !   IF (td(blki(1, blk)) .EQ. ithread) my_cnt = my_cnt+1
         !ENDDO
         thr_c(ithread + 1) = my_cnt
!$OMP        BARRIER
!$OMP        MASTER
         CALL dbcsr_build_row_index_inplace(thr_c, nthreads)
!$OMP        END MASTER
!$OMP        BARRIER
         my_cnt = (thr_c(ithread + 1) + 1)*3 - 2
         IF (matrix%local_indexing) THEN
            DO blk = 1, nblks
               IF (td(local_rows(blki(1, blk))) .EQ. ithread) THEN
                  matrix%coo_l(my_cnt:my_cnt + 2) = blki(1:3, blk)
                  my_cnt = my_cnt + 3
               END IF
            END DO
         ELSE
            DO blk = 1, nblks
               IF (td(blki(1, blk)) .EQ. ithread) THEN
                  matrix%coo_l(my_cnt:my_cnt + 2) = blki(1:3, blk)
                  my_cnt = my_cnt + 3
               END IF
            END DO
         END IF
!$OMP        END PARALLEL
      ELSE
         ! Small price to pay for avoiding infinite recursions.
         DO blk = 2, nblks
            IF (blki(1, blk) .EQ. blki(1, blk - 1) .AND. blki(2, blk) .EQ. blki(2, blk - 1)) THEN
               ! Weird assertion to give some idea of the two blocks.
               IF (-blki(1, blk) /= blki(2, blk)) &
                  CALL dbcsr_abort(__LOCATION__, &
                                   "Should not have duplicate matrix blocks. (-row, col) is duplicated.")
            END IF
         END DO
         !
         IF (nblks > 0) THEN
            CALL dbcsr_addto_index_array(matrix, dbcsr_slot_coo_l, &
                                         DATA=RESHAPE(blki, (/3*nblks/)))
         ELSE
            CALL dbcsr_addto_index_array(matrix, dbcsr_slot_coo_l, &
                                         DATA=zero_len_array)
         END IF
      END IF

      DEALLOCATE (blki)
      !
      CALL timestop(error_handle)
   END SUBROUTINE dbcsr_make_index_list