her2k

Performs a hermitian rank-2k update.

Description

The her2k routines perform a rank-2k update of an n x n hermitian matrix C by general matrices A and B. The operation is defined as:

If trans = transpose::nontrans,

C \leftarrow alpha*A*B^H + conjg(alpha)B*A^H + beta*C

where:

  • A is n x k and B is k x n.

If trans = transpose::conjtrans,

C \leftarrow alpha*A^H*B + conjg(alpha)*B^H*A + beta*C

where:

  • A is k x n and B is n x k.

In both cases:

  • alpha is a complex scalar and beta is a real scalar

  • C is hermitian matrix, A and B are general matrices

  • The inner dimension of both matrix multiplications is k

her2k supports the following precisions:

T

T_real

std::complex<float>

float

std::complex<double>

double

her2k (Buffer Version)

Syntax

namespace oneapi::mkl::blas::column_major {
    void her2k(sycl::queue &queue,
               oneapi::mkl::uplo upper_lower,
               oneapi::mkl::transpose trans,
               std::int64_t n,
               std::int64_t k,
               T alpha,
               sycl::buffer<T,1> &a,
               std::int64_t lda,
               sycl::buffer<T,1> &b,
               std::int64_t ldb,
               T_real beta,
               sycl::buffer<T,1> &c,
               std::int64_t ldc)
}
namespace oneapi::mkl::blas::row_major {
    void her2k(sycl::queue &queue,
               oneapi::mkl::uplo upper_lower,
               oneapi::mkl::transpose trans,
               std::int64_t n,
               std::int64_t k,
               T alpha,
               sycl::buffer<T,1> &a,
               std::int64_t lda,
               sycl::buffer<T,1> &b,
               std::int64_t ldb,
               T_real beta,
               sycl::buffer<T,1> &c,
               std::int64_t ldc)
}

Input Parameters

queue

The queue where the routine should be executed.

upper_lower

Specifies whether matrix C is upper or lower triangular. See Data Types for more details.

trans

Specifies the transposition operation applied as described above. Supported operations are transpose::nontrans and transpose::conjtrans.

n

Number of rows and columns of matrix C. Must be at least zero.

k

Inner dimension of matrix multiplications. Must be at least zero.

alpha

Complex scaling factor for the rank-2k update.

a

Buffer holding input matrix A. See Matrix Storage for more details.

trans = transpose::nontrans

trans = transpose::conjtrans

Column major

A is n x k matrix. Size of array a must be at least lda * k

A is k x n matrix. Size of array a must be at least lda * n

Row major

A is n x k matrix. Size of array a must be at least lda * n

A is k x n matrix. Size of array a must be at least lda * k

lda

Leading dimension of matrix A. Must be positive.

trans = transpose::nontrans

trans = transpose::conjtrans

Column major

Must be at least n

Must be at least k

Row major

Must be at least k

Must be at least n

beta

Real scaling factor for matrix C.

b

Buffer holding input matrix B. See Matrix Storage for more details.

trans = transpose::nontrans

trans = transpose::conjtrans

Column major

B is k x n matrix. Size of array b must be at least ldb * n

B is n x k matrix. Size of array b must be at least ldb * k

Row major

B is k x n matrix. Size of array b must be at least ldb * k

B is n x k matrix. Size of array b must be at least ldb * n

ldb

Leading dimension of matrix B. Must be positive.

trans = transpose::nontrans

trans = transpose::conjtrans

Column major

Must be at least k

Must be at least n

Row major

Must be at least n

Must be at least k

c

Buffer holding input/output matrix C. Size of the buffer must be at least ldc * n. See Matrix Storage for more details.

ldc

Leading dimension of matrix C. Must be positive and at least n.

Output Parameters

c

Output buffer overwritten by updated C matrix.

her2k (USM Version)

Syntax

namespace oneapi::mkl::blas::column_major {
    sycl::event her2k(sycl::queue &queue,
                      oneapi::mkl::uplo upper_lower,
                      oneapi::mkl::transpose trans,
                      std::int64_t n,
                      std::int64_t k,
                      T alpha,
                      const T* a,
                      std::int64_t lda,
                      const T* b,
                      std::int64_t ldb,
                      T_real beta,
                      T* c,
                      std::int64_t ldc,
                      const std::vector<sycl::event> &dependencies = {})
}
namespace oneapi::mkl::blas::row_major {
    sycl::event her2k(sycl::queue &queue,
                      oneapi::mkl::uplo upper_lower,
                      oneapi::mkl::transpose trans,
                      std::int64_t n,
                      std::int64_t k,
                      T alpha,
                      const T* a,
                      std::int64_t lda,
                      const T* b,
                      std::int64_t ldb,
                      T_real beta,
                      T* c,
                      std::int64_t ldc,
                      const std::vector<sycl::event> &dependencies = {})
}

Input Parameters

queue

The queue where the routine should be executed.

upper_lower

Specifies whether matrix C is upper or lower triangular. See Data Types for more details.

trans

Specifies the transposition operation applied as described above. Supported operations are transpose::nontrans and transpose::conjtrans.

n

Number of rows and columns of matrix C. Must be at least zero.

k

Inner dimension of matrix multiplications. Must be at least zero.

alpha

Complex scaling factor for the rank-2k update.

a

Pointer to input matrix A. See Matrix Storage for more details.

trans = transpose::nontrans

trans = transpose::conjtrans

Column major

A is n x k matrix. Size of array a must be at least lda * k

A is k x n matrix. Size of array a must be at least lda * n

Row major

A is n x k matrix. Size of array a must be at least lda * n

A is k x n matrix. Size of array a must be at least lda * k

lda

Leading dimension of matrix A. Must be positive.

trans = transpose::nontrans

trans = transpose::conjtrans

Column major

Must be at least n

Must be at least k

Row major

Must be at least k

Must be at least n

beta

Real scaling factor for matrix C.

b

Pointer to input matrix B. See Matrix Storage for more details.

trans = transpose::nontrans

trans = transpose::conjtrans

Column major

B is k x n matrix. Size of array b must be at least ldb * n

B is n x k matrix. Size of array b must be at least ldb * k

Row major

B is k x n matrix. Size of array b must be at least ldb * k

B is n x k matrix. Size of array b must be at least ldb * n

ldb

Leading dimension of matrix B. Must be positive.

trans = transpose::nontrans

trans = transpose::conjtrans

Column major

Must be at least k

Must be at least n

Row major

Must be at least n

Must be at least k

c

Pointer to input/output matrix C. Size of the array must be at least ldc * n. See Matrix Storage for more details.

ldc

Leading dimension of matrix C. Must be positive and at least n.

dependencies

List of events to wait for before starting computation, if any. If omitted, defaults to no dependencies.

Output Parameters

c

Pointer to output matrix overwritten by updated C matrix.

Return Values

Output event to wait on to ensure computation is complete.