mirror of
https://github.com/cp2k/cp2k.git
synced 2026-07-28 14:15:19 -04:00
Fixes and out-of-place trsm/trmm
- Removed cp_fm_cholesky_restore from cp_fm_basic_linalg - Thereby reduced duplicated code (cp_fm_cholesky_restore). - Prefer cp_fm_cholesky to provide this routine. - Fixed cp_fm_cholesky_restore for pos=right. - Fixed cp_fm_triangular_multiply for side=R. - Removed unused variables. - Adjusted some INTENTs.
This commit is contained in:
Hans Pabst
parent
fde9df4dde
commit
c08dc9e0b6
4 changed files with 93 additions and 189 deletions
|
|
@ -76,3 +76,10 @@
|
|||
#else
|
||||
#define OMP_DEFAULT_NONE_WITH_OOP NONE
|
||||
#endif
|
||||
|
||||
#if defined (__MKL)
|
||||
! Preprocessing is enabled by default, and below header is not language specific
|
||||
! Defines __INTEL_MKL__ (2025), __INTEL_MKL_MINOR__ (0), __INTEL_MKL_UPDATE__ (2),
|
||||
! and __INTEL_MKL_PATCH__ (0) as well INTEL_MKL_VERSION, e.g., 20250200.
|
||||
#include <mkl_version.h>
|
||||
#endif
|
||||
|
|
|
|||
|
|
@ -20,10 +20,10 @@ MODULE cp_dbcsr_cholesky
|
|||
dbcsr_type
|
||||
USE cp_dbcsr_operations, ONLY: copy_dbcsr_to_fm,&
|
||||
copy_fm_to_dbcsr
|
||||
USE cp_fm_basic_linalg, ONLY: cp_fm_cholesky_restore,&
|
||||
cp_fm_potrf,&
|
||||
USE cp_fm_basic_linalg, ONLY: cp_fm_potrf,&
|
||||
cp_fm_potri,&
|
||||
cp_fm_uplo_to_full
|
||||
USE cp_fm_cholesky, ONLY: cp_fm_cholesky_restore
|
||||
USE cp_fm_struct, ONLY: cp_fm_struct_create,&
|
||||
cp_fm_struct_release,&
|
||||
cp_fm_struct_type
|
||||
|
|
|
|||
|
|
@ -64,7 +64,6 @@ MODULE cp_fm_basic_linalg
|
|||
cp_fm_potri, & ! Invert triangular matrix
|
||||
cp_fm_rot_rows, & ! rotates two rows
|
||||
cp_fm_rot_cols, & ! rotates two columns
|
||||
cp_fm_cholesky_restore, & ! apply Cholesky decomposition
|
||||
cp_fm_Gram_Schmidt_orthonorm, & ! Gram-Schmidt orthonormalization of columns of a full matrix, &
|
||||
cp_fm_det, & ! determinant of a real matrix with correct sign
|
||||
cp_fm_matvec ! matrix-vector multiplication (vector replicated)
|
||||
|
|
@ -716,7 +715,7 @@ CONTAINS
|
|||
CHARACTER(len=*), PARAMETER :: routineN = 'cp_fm_schur_product'
|
||||
|
||||
INTEGER :: handle, icol_local, irow_local, mypcol, &
|
||||
myprow, ncol_local, npcol, nprow, &
|
||||
myprow, ncol_local, &
|
||||
nrow_local
|
||||
REAL(KIND=dp), DIMENSION(:, :), POINTER :: a, b, c
|
||||
TYPE(cp_blacs_env_type), POINTER :: context
|
||||
|
|
@ -726,8 +725,6 @@ CONTAINS
|
|||
context => matrix_a%matrix_struct%context
|
||||
myprow = context%mepos(1)
|
||||
mypcol = context%mepos(2)
|
||||
nprow = context%num_pe(1)
|
||||
npcol = context%num_pe(2)
|
||||
|
||||
a => matrix_a%local_data
|
||||
b => matrix_b%local_data
|
||||
|
|
@ -767,7 +764,7 @@ CONTAINS
|
|||
CHARACTER(len=*), PARAMETER :: routineN = 'cp_fm_trace_a0b0t0'
|
||||
|
||||
INTEGER :: handle, mypcol, myprow, ncol_local, &
|
||||
npcol, nprow, nrow_local
|
||||
nrow_local
|
||||
REAL(KIND=dp), DIMENSION(:, :), POINTER :: a, b
|
||||
REAL(KIND=sp), DIMENSION(:, :), POINTER :: a_sp, b_sp
|
||||
TYPE(cp_blacs_env_type), POINTER :: context
|
||||
|
|
@ -778,8 +775,6 @@ CONTAINS
|
|||
context => matrix_a%matrix_struct%context
|
||||
myprow = context%mepos(1)
|
||||
mypcol = context%mepos(2)
|
||||
nprow = context%num_pe(1)
|
||||
npcol = context%num_pe(2)
|
||||
|
||||
group = matrix_a%matrix_struct%para_env
|
||||
|
||||
|
|
@ -1117,7 +1112,7 @@ CONTAINS
|
|||
CHARACTER(len=*), PARAMETER :: routineN = 'cp_fm_triangular_multiply'
|
||||
|
||||
CHARACTER :: side_char, transa, unit_diag, uplo
|
||||
INTEGER :: handle, m, n
|
||||
INTEGER :: handle, mdim, m, n
|
||||
LOGICAL :: invert
|
||||
REAL(KIND=dp) :: al
|
||||
|
||||
|
|
@ -1130,6 +1125,7 @@ CONTAINS
|
|||
al = 1.0_dp
|
||||
CALL cp_fm_get_info(matrix_b, nrow_global=m, ncol_global=n)
|
||||
IF (PRESENT(side)) side_char = side
|
||||
mdim = MERGE(1, 2, 'L' .EQ. side_char)
|
||||
IF (PRESENT(invert_tr)) invert = invert_tr
|
||||
IF (PRESENT(uplo_tr)) uplo = uplo_tr
|
||||
IF (PRESENT(unit_diag_tr)) THEN
|
||||
|
|
@ -1161,7 +1157,7 @@ CONTAINS
|
|||
#else
|
||||
CALL dtrsm(side_char, uplo, transa, unit_diag, m, n, al, &
|
||||
triangular_matrix%local_data(1, 1), &
|
||||
SIZE(triangular_matrix%local_data, 1), &
|
||||
SIZE(triangular_matrix%local_data, mdim), &
|
||||
matrix_b%local_data(1, 1), SIZE(matrix_b%local_data, 1))
|
||||
#endif
|
||||
|
||||
|
|
@ -1176,7 +1172,7 @@ CONTAINS
|
|||
#else
|
||||
CALL dtrmm(side_char, uplo, transa, unit_diag, m, n, al, &
|
||||
triangular_matrix%local_data(1, 1), &
|
||||
SIZE(triangular_matrix%local_data, 1), &
|
||||
SIZE(triangular_matrix%local_data, mdim), &
|
||||
matrix_b%local_data(1, 1), SIZE(matrix_b%local_data, 1))
|
||||
#endif
|
||||
|
||||
|
|
@ -1285,7 +1281,7 @@ CONTAINS
|
|||
CHARACTER :: myuplo
|
||||
INTEGER :: handle, icol_global, irow_global, &
|
||||
mypcol, myprow, ncol_global, &
|
||||
npcol, nprow, nrow_global
|
||||
nrow_global
|
||||
REAL(KIND=dp), DIMENSION(:, :), POINTER :: a
|
||||
REAL(KIND=sp), DIMENSION(:, :), POINTER :: a_sp
|
||||
TYPE(cp_blacs_env_type), POINTER :: context
|
||||
|
|
@ -1315,8 +1311,6 @@ CONTAINS
|
|||
context => matrix%matrix_struct%context
|
||||
myprow = context%mepos(1)
|
||||
mypcol = context%mepos(2)
|
||||
nprow = context%num_pe(1)
|
||||
npcol = context%num_pe(2)
|
||||
|
||||
#if defined(__parallel)
|
||||
|
||||
|
|
@ -2711,96 +2705,6 @@ CONTAINS
|
|||
CPASSERT(info == 0)
|
||||
END SUBROUTINE cp_fm_potri
|
||||
|
||||
! **************************************************************************************************
|
||||
!> \brief ...
|
||||
!> \param fm_matrix ...
|
||||
!> \param neig ...
|
||||
!> \param fm_matrixb ...
|
||||
!> \param fm_matrixout ...
|
||||
!> \param op ...
|
||||
!> \param pos ...
|
||||
!> \param transa ...
|
||||
! **************************************************************************************************
|
||||
SUBROUTINE cp_fm_cholesky_restore(fm_matrix, neig, fm_matrixb, fm_matrixout, op, pos, transa)
|
||||
TYPE(cp_fm_type) :: fm_matrix
|
||||
TYPE(cp_fm_type) :: fm_matrixb
|
||||
TYPE(cp_fm_type) :: fm_matrixout
|
||||
INTEGER, INTENT(IN) :: neig
|
||||
CHARACTER(LEN=*), INTENT(IN) :: op
|
||||
CHARACTER(LEN=*), INTENT(IN) :: pos
|
||||
CHARACTER(LEN=*), INTENT(IN) :: transa
|
||||
|
||||
REAL(KIND=dp), DIMENSION(:, :), POINTER :: a, b, outm
|
||||
REAL(KIND=sp), DIMENSION(:, :), POINTER :: a_sp, b_sp, outm_sp
|
||||
INTEGER :: n, itype
|
||||
REAL(KIND=dp) :: alpha
|
||||
#if defined(__parallel)
|
||||
INTEGER :: i
|
||||
INTEGER, DIMENSION(9) :: desca, descb, descout
|
||||
#endif
|
||||
|
||||
! notice b is the cholesky guy
|
||||
a => fm_matrix%local_data
|
||||
b => fm_matrixb%local_data
|
||||
outm => fm_matrixout%local_data
|
||||
a_sp => fm_matrix%local_data_sp
|
||||
b_sp => fm_matrixb%local_data_sp
|
||||
outm_sp => fm_matrixout%local_data_sp
|
||||
|
||||
n = fm_matrix%matrix_struct%nrow_global
|
||||
itype = 1
|
||||
|
||||
#if defined(__parallel)
|
||||
desca(:) = fm_matrix%matrix_struct%descriptor(:)
|
||||
descb(:) = fm_matrixb%matrix_struct%descriptor(:)
|
||||
descout(:) = fm_matrixout%matrix_struct%descriptor(:)
|
||||
alpha = 1.0_dp
|
||||
DO i = 1, neig
|
||||
IF (fm_matrix%use_sp) THEN
|
||||
CALL pscopy(n, a_sp(1, 1), 1, i, desca, 1, outm_sp(1, 1), 1, i, descout, 1)
|
||||
ELSE
|
||||
CALL pdcopy(n, a(1, 1), 1, i, desca, 1, outm(1, 1), 1, i, descout, 1)
|
||||
END IF
|
||||
END DO
|
||||
IF (op .EQ. "SOLVE") THEN
|
||||
IF (fm_matrix%use_sp) THEN
|
||||
CALL pstrsm(pos, 'U', transa, 'N', n, neig, REAL(alpha, sp), b_sp(1, 1), 1, 1, descb, &
|
||||
outm_sp(1, 1), 1, 1, descout)
|
||||
ELSE
|
||||
CALL pdtrsm(pos, 'U', transa, 'N', n, neig, alpha, b(1, 1), 1, 1, descb, outm(1, 1), 1, 1, descout)
|
||||
END IF
|
||||
ELSE
|
||||
IF (fm_matrix%use_sp) THEN
|
||||
CALL pstrmm(pos, 'U', transa, 'N', n, neig, REAL(alpha, sp), b_sp(1, 1), 1, 1, descb, &
|
||||
outm_sp(1, 1), 1, 1, descout)
|
||||
ELSE
|
||||
CALL pdtrmm(pos, 'U', transa, 'N', n, neig, alpha, b(1, 1), 1, 1, descb, outm(1, 1), 1, 1, descout)
|
||||
END IF
|
||||
END IF
|
||||
#else
|
||||
alpha = 1.0_dp
|
||||
IF (fm_matrix%use_sp) THEN
|
||||
CALL scopy(neig*n, a_sp(1, 1), 1, outm_sp(1, 1), 1)
|
||||
ELSE
|
||||
CALL dcopy(neig*n, a(1, 1), 1, outm(1, 1), 1)
|
||||
END IF
|
||||
IF (op .EQ. "SOLVE") THEN
|
||||
IF (fm_matrix%use_sp) THEN
|
||||
CALL strsm(pos, 'U', transa, 'N', n, neig, REAL(alpha, sp), b_sp(1, 1), SIZE(b_sp, 1), outm_sp(1, 1), n)
|
||||
ELSE
|
||||
CALL dtrsm(pos, 'U', transa, 'N', n, neig, alpha, b(1, 1), SIZE(b, 1), outm(1, 1), n)
|
||||
END IF
|
||||
ELSE
|
||||
IF (fm_matrix%use_sp) THEN
|
||||
CALL strmm(pos, 'U', transa, 'N', n, neig, REAL(alpha, sp), b_sp(1, 1), n, outm_sp(1, 1), n)
|
||||
ELSE
|
||||
CALL dtrmm(pos, 'U', transa, 'N', n, neig, alpha, b(1, 1), n, outm(1, 1), n)
|
||||
END IF
|
||||
END IF
|
||||
#endif
|
||||
|
||||
END SUBROUTINE cp_fm_cholesky_restore
|
||||
|
||||
! **************************************************************************************************
|
||||
!> \brief Calculates
|
||||
!> yv = alpha*amat*xv + beta*yv
|
||||
|
|
|
|||
|
|
@ -12,7 +12,6 @@
|
|||
!> \author Fawzi Mohamed
|
||||
! **************************************************************************************************
|
||||
MODULE cp_fm_cholesky
|
||||
USE cp_blacs_env, ONLY: cp_blacs_env_type
|
||||
USE cp_dlaf_utils_api, ONLY: cp_dlaf_create_grid,&
|
||||
cp_dlaf_initialize
|
||||
USE cp_fm_dlaf_api, ONLY: cp_pdpotrf_dlaf,&
|
||||
|
|
@ -140,13 +139,13 @@ CONTAINS
|
|||
!> \author Joost VandeVondele
|
||||
! **************************************************************************************************
|
||||
SUBROUTINE cp_fm_cholesky_invert(matrix, n, info_out)
|
||||
TYPE(cp_fm_type), INTENT(IN) :: matrix
|
||||
TYPE(cp_fm_type), INTENT(IN) :: matrix
|
||||
INTEGER, INTENT(in), OPTIONAL :: n
|
||||
INTEGER, INTENT(OUT), OPTIONAL :: info_out
|
||||
|
||||
CHARACTER(len=*), PARAMETER :: routineN = 'cp_fm_cholesky_invert'
|
||||
REAL(KIND=dp), DIMENSION(:, :), POINTER :: a
|
||||
REAL(KIND=sp), DIMENSION(:, :), POINTER :: a_sp
|
||||
REAL(KIND=dp), DIMENSION(:, :), POINTER :: a
|
||||
REAL(KIND=sp), DIMENSION(:, :), POINTER :: a_sp
|
||||
INTEGER :: info, handle
|
||||
INTEGER :: my_n
|
||||
#if defined(__parallel)
|
||||
|
|
@ -228,11 +227,11 @@ CONTAINS
|
|||
!> \author Joost VandeVondele
|
||||
! **************************************************************************************************
|
||||
SUBROUTINE cp_fm_cholesky_reduce(matrix, matrixb, itype)
|
||||
TYPE(cp_fm_type), INTENT(IN) :: matrix, matrixb
|
||||
INTEGER, OPTIONAL :: itype
|
||||
TYPE(cp_fm_type), INTENT(IN) :: matrix, matrixb
|
||||
INTEGER, OPTIONAL :: itype
|
||||
|
||||
CHARACTER(len=*), PARAMETER :: routineN = 'cp_fm_cholesky_reduce'
|
||||
REAL(KIND=dp), DIMENSION(:, :), POINTER :: a, b
|
||||
REAL(KIND=dp), DIMENSION(:, :), POINTER :: a, b
|
||||
INTEGER :: info, handle
|
||||
INTEGER :: n, my_itype
|
||||
#if defined(__parallel)
|
||||
|
|
@ -273,130 +272,124 @@ CONTAINS
|
|||
|
||||
END SUBROUTINE cp_fm_cholesky_reduce
|
||||
|
||||
!
|
||||
! op can be "SOLVE" (out = U^-1 * in ) or "MULTIPLY" (out = U * in )
|
||||
! pos can be "LEFT" or "RIGHT" (U at the left or at the right)
|
||||
!
|
||||
! DEPRECATED, see cp_fm_basic_linalg:cp_fm_triangular_multiply
|
||||
!
|
||||
! **************************************************************************************************
|
||||
!> \brief ...
|
||||
!> \param matrix ...
|
||||
!> \brief apply Cholesky decomposition
|
||||
!> op can be "SOLVE" (out = U^-1 * in) or "MULTIPLY" (out = U * in)
|
||||
!> pos can be "LEFT" or "RIGHT" (U at the left or at the right)
|
||||
!> \param fm_matrix ...
|
||||
!> \param neig ...
|
||||
!> \param matrixb ...
|
||||
!> \param matrixout ...
|
||||
!> \param fm_matrixb ...
|
||||
!> \param fm_matrixout ...
|
||||
!> \param op ...
|
||||
!> \param pos ...
|
||||
!> \param transa ...
|
||||
! **************************************************************************************************
|
||||
SUBROUTINE cp_fm_cholesky_restore(matrix, neig, matrixb, matrixout, op, pos, transa)
|
||||
TYPE(cp_fm_type), INTENT(IN) :: matrix, matrixb, matrixout
|
||||
INTEGER, INTENT(IN) :: neig
|
||||
CHARACTER(LEN=*), INTENT(IN) :: op
|
||||
CHARACTER(LEN=*), INTENT(IN), OPTIONAL :: pos
|
||||
CHARACTER(LEN=*), INTENT(IN), OPTIONAL :: transa
|
||||
SUBROUTINE cp_fm_cholesky_restore(fm_matrix, neig, fm_matrixb, fm_matrixout, op, pos, transa)
|
||||
TYPE(cp_fm_type), INTENT(IN) :: fm_matrix, fm_matrixb, fm_matrixout
|
||||
INTEGER, INTENT(IN) :: neig
|
||||
CHARACTER(LEN=*), INTENT(IN) :: op
|
||||
CHARACTER(LEN=*), INTENT(IN), OPTIONAL :: pos
|
||||
CHARACTER(LEN=*), INTENT(IN), OPTIONAL :: transa
|
||||
|
||||
CHARACTER(len=*), PARAMETER :: routineN = 'cp_fm_cholesky_restore'
|
||||
REAL(KIND=dp), DIMENSION(:, :), POINTER :: a, b, out
|
||||
REAL(KIND=sp), DIMENSION(:, :), POINTER :: a_sp, b_sp, out_sp
|
||||
INTEGER :: itype, handle
|
||||
INTEGER :: n
|
||||
REAL(KIND=dp) :: alpha
|
||||
INTEGER :: myprow, mypcol
|
||||
TYPE(cp_blacs_env_type), POINTER :: context
|
||||
REAL(KIND=dp), DIMENSION(:, :), POINTER :: a, b, outm
|
||||
REAL(KIND=sp), DIMENSION(:, :), POINTER :: a_sp, b_sp, outm_sp
|
||||
REAL(KIND=dp) :: alpha
|
||||
INTEGER :: handle, mdim, n
|
||||
CHARACTER :: chol_pos, chol_transa
|
||||
#if defined(__parallel)
|
||||
INTEGER :: i
|
||||
INTEGER, DIMENSION(9) :: desca, descb, descout
|
||||
#elif defined(__MKL) && (2025 <= __INTEL_MKL__)
|
||||
REAL(KIND=dp) :: beta
|
||||
#endif
|
||||
|
||||
CALL timeset(routineN, handle)
|
||||
|
||||
context => matrix%matrix_struct%context
|
||||
myprow = context%mepos(1)
|
||||
mypcol = context%mepos(2)
|
||||
n = matrix%matrix_struct%nrow_global
|
||||
itype = 1
|
||||
IF (op /= "SOLVE" .AND. op /= "MULTIPLY") &
|
||||
CPABORT("wrong argument op")
|
||||
! wrong argument op
|
||||
CPASSERT((op == "SOLVE") .OR. (op == "MULTIPLY"))
|
||||
! not the same precision
|
||||
CPASSERT(fm_matrix%use_sp .EQV. fm_matrixb%use_sp)
|
||||
CPASSERT(fm_matrix%use_sp .EQV. fm_matrixout%use_sp)
|
||||
|
||||
IF (PRESENT(pos)) THEN
|
||||
SELECT CASE (pos)
|
||||
CASE ("LEFT")
|
||||
chol_pos = 'L'
|
||||
CASE ("RIGHT")
|
||||
chol_pos = 'R'
|
||||
CASE DEFAULT
|
||||
CPABORT("wrong argument pos")
|
||||
END SELECT
|
||||
ELSE
|
||||
chol_pos = 'L'
|
||||
END IF
|
||||
chol_pos = 'L'
|
||||
IF (PRESENT(pos)) chol_pos = pos(1:1)
|
||||
CPASSERT((chol_pos == 'L') .OR. (chol_pos == 'R'))
|
||||
|
||||
chol_transa = 'N'
|
||||
IF (PRESENT(transa)) chol_transa = transa
|
||||
|
||||
IF ((matrix%use_sp .NEQV. matrixb%use_sp) .OR. (matrix%use_sp .NEQV. matrixout%use_sp)) &
|
||||
CPABORT("not the same precision")
|
||||
|
||||
! notice b is the cholesky guy
|
||||
a => matrix%local_data
|
||||
b => matrixb%local_data
|
||||
out => matrixout%local_data
|
||||
a_sp => matrix%local_data_sp
|
||||
b_sp => matrixb%local_data_sp
|
||||
out_sp => matrixout%local_data_sp
|
||||
a => fm_matrix%local_data
|
||||
b => fm_matrixb%local_data
|
||||
outm => fm_matrixout%local_data
|
||||
a_sp => fm_matrix%local_data_sp
|
||||
b_sp => fm_matrixb%local_data_sp
|
||||
outm_sp => fm_matrixout%local_data_sp
|
||||
n = fm_matrix%matrix_struct%nrow_global
|
||||
mdim = MERGE(1, 2, 'L' .EQ. chol_pos)
|
||||
alpha = 1.0_dp
|
||||
|
||||
#if defined(__parallel)
|
||||
|
||||
desca(:) = matrix%matrix_struct%descriptor(:)
|
||||
descb(:) = matrixb%matrix_struct%descriptor(:)
|
||||
descout(:) = matrixout%matrix_struct%descriptor(:)
|
||||
alpha = 1.0_dp
|
||||
desca(:) = fm_matrix%matrix_struct%descriptor(:)
|
||||
descb(:) = fm_matrixb%matrix_struct%descriptor(:)
|
||||
descout(:) = fm_matrixout%matrix_struct%descriptor(:)
|
||||
DO i = 1, neig
|
||||
IF (matrix%use_sp) THEN
|
||||
CALL pscopy(n, a_sp(1, 1), 1, i, desca, 1, out_sp(1, 1), 1, i, descout, 1)
|
||||
IF (fm_matrix%use_sp) THEN
|
||||
CALL pscopy(n, a_sp(1, 1), 1, i, desca, 1, outm_sp(1, 1), 1, i, descout, 1)
|
||||
ELSE
|
||||
CALL pdcopy(n, a(1, 1), 1, i, desca, 1, out(1, 1), 1, i, descout, 1)
|
||||
CALL pdcopy(n, a(1, 1), 1, i, desca, 1, outm(1, 1), 1, i, descout, 1)
|
||||
END IF
|
||||
END DO
|
||||
IF (op .EQ. "SOLVE") THEN
|
||||
IF (matrix%use_sp) THEN
|
||||
IF (fm_matrix%use_sp) THEN
|
||||
CALL pstrsm(chol_pos, 'U', chol_transa, 'N', n, neig, REAL(alpha, sp), b_sp(1, 1), 1, 1, descb, &
|
||||
out_sp(1, 1), 1, 1, descout)
|
||||
outm_sp(1, 1), 1, 1, descout)
|
||||
ELSE
|
||||
CALL pdtrsm(chol_pos, 'U', chol_transa, 'N', n, neig, alpha, b(1, 1), 1, 1, descb, out(1, 1), 1, 1, descout)
|
||||
CALL pdtrsm(chol_pos, 'U', chol_transa, 'N', n, neig, alpha, b(1, 1), 1, 1, descb, outm(1, 1), 1, 1, descout)
|
||||
END IF
|
||||
ELSE
|
||||
IF (matrix%use_sp) THEN
|
||||
IF (fm_matrix%use_sp) THEN
|
||||
CALL pstrmm(chol_pos, 'U', chol_transa, 'N', n, neig, REAL(alpha, sp), b_sp(1, 1), 1, 1, descb, &
|
||||
out_sp(1, 1), 1, 1, descout)
|
||||
outm_sp(1, 1), 1, 1, descout)
|
||||
ELSE
|
||||
CALL pdtrmm(chol_pos, 'U', chol_transa, 'N', n, neig, alpha, b(1, 1), 1, 1, descb, out(1, 1), 1, 1, descout)
|
||||
CALL pdtrmm(chol_pos, 'U', chol_transa, 'N', n, neig, alpha, b(1, 1), 1, 1, descb, outm(1, 1), 1, 1, descout)
|
||||
END IF
|
||||
END IF
|
||||
#elif defined(__MKL) && (2025 <= __INTEL_MKL__)
|
||||
beta = 0.0_dp
|
||||
IF (op .EQ. "SOLVE") THEN
|
||||
IF (fm_matrix%use_sp) THEN
|
||||
CALL strsm_oop(chol_pos, 'U', chol_transa, 'N', n, neig, REAL(alpha, sp), b_sp(1, 1), SIZE(b_sp, mdim), a_sp(1, 1), n, REAL(beta, sp), outm_sp(1, 1), n)
|
||||
ELSE
|
||||
CALL dtrsm_oop(chol_pos, 'U', chol_transa, 'N', n, neig, alpha, b(1, 1), SIZE(b, mdim), a(1, 1), n, beta, outm(1, 1), n)
|
||||
END IF
|
||||
ELSE
|
||||
IF (fm_matrix%use_sp) THEN
|
||||
CALL strmm_oop(chol_pos, 'U', chol_transa, 'N', n, neig, REAL(alpha, sp), b_sp(1, 1), SIZE(b_sp, mdim), a_sp(1, 1), n, REAL(beta, sp), outm_sp(1, 1), n)
|
||||
ELSE
|
||||
CALL dtrmm_oop(chol_pos, 'U', chol_transa, 'N', n, neig, alpha, b(1, 1), SIZE(b, mdim), a(1, 1), n, beta, outm(1, 1), n)
|
||||
END IF
|
||||
END IF
|
||||
#else
|
||||
|
||||
alpha = 1.0_dp
|
||||
IF (matrix%use_sp) THEN
|
||||
CALL scopy(neig*n, a_sp(1, 1), 1, out_sp(1, 1), 1)
|
||||
IF (fm_matrix%use_sp) THEN
|
||||
CALL scopy(neig*n, a_sp(1, 1), 1, outm_sp(1, 1), 1)
|
||||
ELSE
|
||||
CALL dcopy(neig*n, a(1, 1), 1, out(1, 1), 1)
|
||||
CALL dcopy(neig*n, a(1, 1), 1, outm(1, 1), 1)
|
||||
END IF
|
||||
IF (op .EQ. "SOLVE") THEN
|
||||
IF (matrix%use_sp) THEN
|
||||
CALL strsm(chol_pos, 'U', chol_transa, 'N', n, neig, REAL(alpha, sp), b_sp(1, 1), SIZE(b_sp, 1), out_sp(1, 1), n)
|
||||
IF (fm_matrix%use_sp) THEN
|
||||
CALL strsm(chol_pos, 'U', chol_transa, 'N', n, neig, REAL(alpha, sp), b_sp(1, 1), SIZE(b_sp, mdim), outm_sp(1, 1), n)
|
||||
ELSE
|
||||
CALL dtrsm(chol_pos, 'U', chol_transa, 'N', n, neig, alpha, b(1, 1), SIZE(b, 1), out(1, 1), n)
|
||||
CALL dtrsm(chol_pos, 'U', chol_transa, 'N', n, neig, alpha, b(1, 1), SIZE(b, mdim), outm(1, 1), n)
|
||||
END IF
|
||||
ELSE
|
||||
IF (matrix%use_sp) THEN
|
||||
CALL strmm(chol_pos, 'U', chol_transa, 'N', n, neig, REAL(alpha, sp), b_sp(1, 1), n, out_sp(1, 1), n)
|
||||
IF (fm_matrix%use_sp) THEN
|
||||
CALL strmm(chol_pos, 'U', chol_transa, 'N', n, neig, REAL(alpha, sp), b_sp(1, 1), SIZE(b_sp, mdim), outm_sp(1, 1), n)
|
||||
ELSE
|
||||
CALL dtrmm(chol_pos, 'U', chol_transa, 'N', n, neig, alpha, b(1, 1), n, out(1, 1), n)
|
||||
CALL dtrmm(chol_pos, 'U', chol_transa, 'N', n, neig, alpha, b(1, 1), SIZE(b, mdim), outm(1, 1), n)
|
||||
END IF
|
||||
END IF
|
||||
|
||||
#endif
|
||||
|
||||
CALL timestop(handle)
|
||||
|
|
|
|||
Loading…
Add table
Add a link
Reference in a new issue