dbcsr_make_random_matrix Subroutine

public subroutine dbcsr_make_random_matrix(matrix, row_blk_sizes, col_blk_sizes, name, sparsity, mp_group, data_type, symmetry, dist)

Creates a random matrix.

Arguments

Type IntentOptional Attributes Name
type(dbcsr_type), intent(out) :: matrix
integer, intent(inout), DIMENSION(:), POINTER, CONTIGUOUS :: row_blk_sizes
integer, intent(inout), DIMENSION(:), POINTER, CONTIGUOUS :: col_blk_sizes
character(len=*), intent(in) :: name
real(kind=real_8), intent(in) :: sparsity
type(mp_comm_type), intent(in) :: mp_group
integer, intent(in), optional :: data_type
character(len=1), intent(in), optional :: symmetry
type(dbcsr_distribution_obj), intent(in), optional :: dist

Source Code

   SUBROUTINE dbcsr_make_random_matrix(matrix, row_blk_sizes, col_blk_sizes, &
                                       name, sparsity, mp_group, data_type, symmetry, dist)
      !! Creates a random matrix.
      TYPE(dbcsr_type), INTENT(out)                      :: matrix
      INTEGER, DIMENSION(:), INTENT(INOUT), POINTER, CONTIGUOUS :: row_blk_sizes, col_blk_sizes
      CHARACTER(len=*), INTENT(in)                       :: name
      REAL(kind=real_8), INTENT(in)                      :: sparsity
      TYPE(mp_comm_type), INTENT(in)                     :: mp_group
      INTEGER, INTENT(in), OPTIONAL                      :: data_type
      CHARACTER, INTENT(in), OPTIONAL                    :: symmetry
      TYPE(dbcsr_distribution_obj), INTENT(IN), OPTIONAL :: dist

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

      CHARACTER                                          :: my_symmetry
      INTEGER                                            :: col, error_handle, max_nze, &
                                                            my_data_type, my_proc, ncol, nrow, &
                                                            numproc, nze, p, row, s_col, s_row
      INTEGER(KIND=int_8)                                :: counter, ele, increment, nmax
      INTEGER, DIMENSION(4)                              :: iseed, jseed
      LOGICAL                                            :: tr
      REAL(kind=real_8)                                  :: my_sparsity
      REAL(kind=real_8), DIMENSION(1)                    :: value
      TYPE(dbcsr_data_obj)                               :: data_values, data_values_tr
      TYPE(dbcsr_distribution_obj)                       :: new_dist

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

      CALL timeset(routineN, error_handle)
      ! Check that the counter was initialised (or has not overflowed)
      DBCSR_ASSERT(randmat_counter .NE. 0)
      ! the counter goes into the seed. Every new call gives a new random matrix
      randmat_counter = randmat_counter + 1
      ! Create the matrix
      IF (PRESENT(dist)) THEN
         new_dist = dist
         CALL dbcsr_distribution_hold(new_dist)
      ELSE
         CALL dbcsr_make_null_dist(new_dist, SIZE(row_blk_sizes), &
                                   SIZE(col_blk_sizes), group=mp_group)
      END IF
      my_data_type = dbcsr_type_real_default
      IF (PRESENT(data_type)) my_data_type = data_type
      my_symmetry = dbcsr_type_no_symmetry
      IF (PRESENT(symmetry)) my_symmetry = symmetry
      CALL dbcsr_create(matrix, name, &
                        new_dist, my_symmetry, &
                        row_blk_sizes, &
                        col_blk_sizes, &
                        data_type=my_data_type)
      numproc = dbcsr_mp_numnodes(dbcsr_distribution_mp(new_dist))
      my_proc = dbcsr_mp_mynode(dbcsr_distribution_mp(new_dist))
      !
      IF (sparsity .GT. 1) THEN
         my_sparsity = sparsity/100.0
      ELSE
         my_sparsity = sparsity
      END IF
      CALL dbcsr_work_create(matrix, &
                             nblks_guess=INT(REAL(dbcsr_nblkrows_total(matrix), KIND=dp) &
                                             *REAL(dbcsr_nblkcols_total(matrix), KIND=dp) &
                                             *(1.0_dp - sparsity)*1.1_dp/numproc), &
                             sizedata_guess=INT(REAL(dbcsr_nfullrows_total(matrix), KIND=dp) &
                                                *REAL(dbcsr_nfullcols_total(matrix), KIND=dp) &
                                                *(1.0_dp - sparsity)*1.1_dp/numproc), &
                             work_mutable=.TRUE.)

      max_nze = dbcsr_max_row_size(matrix)*dbcsr_max_col_size(matrix)
      CALL dbcsr_data_init(data_values)
      CALL dbcsr_data_new(data_values, my_data_type, data_size=max_nze)
      CALL dbcsr_data_init(data_values_tr)
      CALL dbcsr_data_new(data_values_tr, my_data_type, data_size=max_nze)

      nrow = dbcsr_nblkrows_total(matrix)
      ncol = dbcsr_nblkcols_total(matrix)
      nmax = INT(nrow, KIND=int_8)*INT(ncol, KIND=int_8)
      ele = -1
      counter = 0
      CALL set_larnv_seed(7, 42, 3, 42, randmat_counter, jseed)

      DO
         ! find the next block to add, this is given by a geometrically distributed variable
         ! we number the blocks of the matrix and jump to the next one
         CALL dlarnv(1, jseed, 1, value)
         IF (my_sparsity > 0) THEN
            increment = 1 + FLOOR(LOG(value(1))/LOG(my_sparsity), KIND=int_8)
         ELSE
            increment = 1
         END IF
         ele = ele + increment
         IF (ele >= nmax) EXIT
         counter = counter + 1
         row = INT(ele/ncol) + 1
         col = INT(MOD(ele, INT(ncol, KIND=KIND(ele)))) + 1

         ! build the upper matrix if some symmetry, and only deal with the local blocks.
         s_row = row; s_col = col
         IF (PRESENT(dist)) THEN
            tr = .FALSE.
            CALL dbcsr_get_stored_coordinates(matrix, s_row, s_col, p)
            IF (my_symmetry .NE. dbcsr_type_no_symmetry .AND. s_col .LT. s_row) CYCLE
            IF (p .NE. my_proc) CYCLE
         ELSE
            IF (my_symmetry .NE. dbcsr_type_no_symmetry .AND. s_col .LT. s_row) CYCLE
         END IF
         IF (.NOT. PRESENT(dist) .AND. my_proc .NE. 0) CYCLE

         ! fill based on a block based seed, makes this the same values in parallel
         CALL set_larnv_seed(row, nrow, col, ncol, randmat_counter, iseed)
         nze = row_blk_sizes(s_row)*col_blk_sizes(s_col)
         CALL dbcsr_lapack_larnv(1, iseed, nze, data_values)
         CALL dbcsr_put_block(matrix, s_row, s_col, data_values)
         IF (my_symmetry .NE. dbcsr_type_no_symmetry .AND. s_col .EQ. s_row) THEN
            SELECT CASE (my_symmetry)
            CASE (dbcsr_type_symmetric)
               CALL dbcsr_block_transpose(data_values_tr, data_values, &
                                          row_size=row_blk_sizes(s_row), col_size=col_blk_sizes(s_col), lb=1, source_lb=1)
            CASE (dbcsr_type_antisymmetric)
               CALL dbcsr_block_transpose(data_values_tr, data_values, &
                                          row_size=row_blk_sizes(s_row), col_size=col_blk_sizes(s_col), lb=1, source_lb=1, &
                                          scale=dbcsr_scalar_negative(dbcsr_scalar_one(my_data_type)))
            CASE (dbcsr_type_hermitian)
               CALL dbcsr_block_transpose(data_values_tr, data_values, &
                                          row_size=row_blk_sizes(s_row), col_size=col_blk_sizes(s_col), lb=1, source_lb=1)
               CALL dbcsr_block_conjg(data_values_tr, row_size=col_blk_sizes(s_col), col_size=row_blk_sizes(s_row), &
                                      lb=1)
            CASE (dbcsr_type_antihermitian)
               CALL dbcsr_block_transpose(data_values_tr, data_values, &
                                          row_size=row_blk_sizes(s_row), col_size=col_blk_sizes(s_col), lb=1, source_lb=1, &
                                          scale=dbcsr_scalar_negative(dbcsr_scalar_one(my_data_type)))
               CALL dbcsr_block_conjg(data_values_tr, row_size=col_blk_sizes(s_col), col_size=row_blk_sizes(s_row), &
                                      lb=1)
            CASE DEFAULT
               DBCSR_ABORT("wrong matrix symmetry")
            END SELECT
            CALL dbcsr_put_block(matrix, s_row, s_col, data_values_tr, summation=.TRUE.)
         END IF
      END DO

      CALL dbcsr_data_release(data_values)
      CALL dbcsr_data_release(data_values_tr)

      CALL dbcsr_distribution_release(new_dist)
      CALL dbcsr_finalize(matrix)
      CALL dbcsr_verify_matrix(matrix)
      !
      CALL timestop(error_handle)
   END SUBROUTINE dbcsr_make_random_matrix