add and scale matrices A = alphaA + betaB or
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
type(dbcsr_type), | intent(inout) | :: | matrix_a |
DBCSR matrix |
||
type(dbcsr_type), | intent(in) | :: | matrix_b |
DBCSR matrix |
||
type(dbcsr_scalar_type), | intent(in), | optional | :: | alpha_scalar | ||
type(dbcsr_scalar_type), | intent(in), | optional | :: | beta_scalar | ||
integer(kind=int_8), | intent(inout), | optional | :: | flop |
SUBROUTINE dbcsr_add_anytype(matrix_a, matrix_b, alpha_scalar, beta_scalar, flop)
!! add and scale matrices
!! A = alpha*A + beta*B or
TYPE(dbcsr_type), INTENT(INOUT) :: matrix_a
!! DBCSR matrix
TYPE(dbcsr_type), INTENT(IN) :: matrix_b
!! DBCSR matrix
TYPE(dbcsr_scalar_type), INTENT(IN), OPTIONAL :: alpha_scalar, beta_scalar
INTEGER(KIND=int_8), INTENT(INOUT), OPTIONAL :: flop
CHARACTER(len=*), PARAMETER :: routineN = 'dbcsr_add_anytype'
INTEGER :: data_type_a, data_type_b, &
handle, size_work, iw
INTEGER(KIND=int_8) :: my_flop, local_matrix_size
LOGICAL :: do_scale
TYPE(dbcsr_iterator) :: iter
TYPE(dbcsr_scalar_type) :: my_beta_scalar
! ---------------------------------------------------------------------------
CALL timeset(routineN, handle)
IF (.NOT. dbcsr_valid_index(matrix_a)) &
DBCSR_ABORT("Invalid matrix")
IF ((dbcsr_get_matrix_type(matrix_a) .EQ. dbcsr_type_symmetric .OR. &
dbcsr_get_matrix_type(matrix_a) .EQ. dbcsr_type_antisymmetric) .NEQV. &
(dbcsr_get_matrix_type(matrix_b) .EQ. dbcsr_type_symmetric .OR. &
dbcsr_get_matrix_type(matrix_b) .EQ. dbcsr_type_antisymmetric)) THEN
DBCSR_ABORT("Summing general with symmetric matrix NYI")
END IF
data_type_a = dbcsr_get_data_type(matrix_a)
data_type_b = dbcsr_get_data_type(matrix_b)
!
my_beta_scalar = dbcsr_scalar_one(data_type_b)
IF (PRESENT(beta_scalar)) my_beta_scalar = beta_scalar
!
! let's go
IF ((dbcsr_nblkrows_total(matrix_a) .NE. dbcsr_nblkrows_total(matrix_b)) .OR. &
(dbcsr_nblkcols_total(matrix_a) .NE. dbcsr_nblkcols_total(matrix_b)) .OR. &
(data_type_a .NE. data_type_b)) &
DBCSR_ABORT("matrices not consistent")
IF (data_type_a .NE. my_beta_scalar%data_type) &
DBCSR_ABORT("beta type parameter not consistent with matrices type")
do_scale = .NOT. dbcsr_scalar_are_equal(my_beta_scalar, dbcsr_scalar_one(data_type_b))
IF (PRESENT(alpha_scalar)) THEN
CALL dbcsr_scale(matrix_a, alpha_scalar=alpha_scalar)
END IF
IF ((.NOT. dbcsr_scalar_are_equal(my_beta_scalar, &
dbcsr_scalar_zero(data_type_b))) .AND. &
dbcsr_get_num_blocks(matrix_b) .GT. 0) THEN
! Pre-size work arrays of matrix_a to avoid continuous reallocation.
! Overestimate for symmetric matrix and multiple threads!
local_matrix_size = INT(dbcsr_nfullrows_local(matrix_a), KIND=int_8)* &
dbcsr_nfullcols_local(matrix_a)
size_work = MAX(0, INT(MIN(local_matrix_size - INT(dbcsr_get_nze(matrix_a), KIND=int_8), &
INT(dbcsr_get_nze(matrix_b), KIND=int_8)), KIND=int_4))
my_flop = 0
!$OMP PARALLEL DEFAULT (NONE) &
!$OMP PRIVATE (iter, iw) &
!$OMP SHARED (matrix_a, matrix_b, data_type_b, size_work) &
!$OMP SHARED (do_scale, my_beta_scalar) &
!$OMP REDUCTION (+ : my_flop)
CALL dbcsr_work_create(matrix_a, &
nblks_guess=matrix_b%nblks, &
sizedata_guess=size_work, &
work_mutable=.FALSE.)
!$OMP BARRIER
iw = 1
!$ iw = omp_get_thread_num() + 1
CALL dbcsr_iterator_start(iter, matrix_b, &
shared=.TRUE., read_only=.TRUE., contiguous_pointers=.FALSE., &
dynamic=.TRUE., dynamic_byrows=.TRUE.)
SELECT CASE (data_type_b)
CASE (dbcsr_type_real_4)
CALL dbcsr_add_anytype_s(matrix_a, matrix_b, iter, iw, do_scale, my_beta_scalar, my_flop)
CASE (dbcsr_type_real_8)
CALL dbcsr_add_anytype_d(matrix_a, matrix_b, iter, iw, do_scale, my_beta_scalar, my_flop)
CASE (dbcsr_type_complex_4)
CALL dbcsr_add_anytype_c(matrix_a, matrix_b, iter, iw, do_scale, my_beta_scalar, my_flop)
CASE (dbcsr_type_complex_8)
CALL dbcsr_add_anytype_z(matrix_a, matrix_b, iter, iw, do_scale, my_beta_scalar, my_flop)
CASE default
DBCSR_ABORT("Invalid data type")
END SELECT
CALL dbcsr_iterator_stop(iter)
CALL dbcsr_finalize(matrix_a)
!$OMP END PARALLEL
IF (PRESENT(flop)) flop = flop + my_flop
END IF
CALL timestop(handle)
END SUBROUTINE dbcsr_add_anytype