.. _onemkl_blas_herk: herk ==== Performs a hermitian rank-k update. Description *********** The ``herk`` routines compute a rank-k update of a hermitian matrix ``C`` by a general matrix ``A``. The operation is defined as: .. math:: C \leftarrow alpha*op(A)*op(A)^H + beta*C where: - op(``X``) is one of op(``X``) = ``X`` or op(``X``) = ``X``\ :sup:`H` - ``alpha`` and ``beta`` are real scalars - ``C`` is ``n`` x ``n`` hermitian matrix - op(``A``) is ``n`` x ``k`` general matrix ``herk`` supports the following precisions: .. list-table:: :header-rows: 1 * - T - T_real * - ``std::complex`` - ``float`` * - ``std::complex`` - ``double`` herk (Buffer Version) ********************* Syntax ------ .. code-block:: cpp namespace oneapi::mkl::blas::column_major { void herk(sycl::queue &queue, oneapi::mkl::uplo upper_lower, oneapi::mkl::transpose trans, std::int64_t n, std::int64_t k, T_real alpha, sycl::buffer &a, std::int64_t lda, T_real beta, sycl::buffer &c, std::int64_t ldc) } .. code-block:: cpp namespace oneapi::mkl::blas::row_major { void herk(sycl::queue &queue, oneapi::mkl::uplo upper_lower, oneapi::mkl::transpose trans, std::int64_t n, std::int64_t k, T_real alpha, sycl::buffer &a, std::int64_t lda, T_real beta, sycl::buffer &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 :ref:`data-types` for more details. trans Specifies op(``A``), the transposition operation applied to matrix ``A``. Supported operations are ``transpose::nontrans`` and ``transpose::conjtrans``. See :ref:`data-types` for more details. n Number of rows and columns of matrix ``C``. Must be at least zero. k Number of columns of matrix op(``A``). Must be at least zero. alpha Complex scaling factor for the rank-k update. a Buffer holding input matrix ``A``. See :ref:`matrix-storage` for more details. .. list-table:: :header-rows: 1 * - - ``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. .. list-table:: :header-rows: 1 * - - ``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``. c Buffer holding input/output matrix ``C``. Size of the buffer must be at least ``ldc`` * ``n``. See :ref:`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 ``alpha`` * op(``A``) * op(``A``)\ :sup:`H` + ``beta`` * ``C``. The imaginary parts of the diagonal elements are set to zero. herk (USM Version) ****************** Syntax ------ .. code-block:: cpp namespace oneapi::mkl::blas::column_major { sycl::event herk(sycl::queue &queue, oneapi::mkl::uplo upper_lower, oneapi::mkl::transpose trans, std::int64_t n, std::int64_t k, T_real alpha, const T* a, std::int64_t lda, T_real beta, T* c, std::int64_t ldc, const std::vector &dependencies = {}) } .. code-block:: cpp namespace oneapi::mkl::blas::row_major { sycl::event herk(sycl::queue &queue, oneapi::mkl::uplo upper_lower, oneapi::mkl::transpose trans, std::int64_t n, std::int64_t k, T_real alpha, const T* a, std::int64_t lda, T_real beta, T* c, std::int64_t ldc, const std::vector &dependencies = {}) } Input Parameters ---------------- queue The queue where the routine should be executed. upper_lower Specifies whether matrix ``C`` is upper or lower triangular. See :ref:`data-types` for more details. trans Specifies op(``A``), the transposition operation applied to matrix ``A``. Supported operations are ``transpose::nontrans`` and ``transpose::conjtrans``. See :ref:`data-types` for more details. n Number of rows and columns of matrix ``C``. Must be at least zero. k Number of columns of matrix op(``A``). Must be at least zero. alpha Complex scaling factor for the rank-k update. a Pointer to input matrix ``A``. See :ref:`matrix-storage` for more details. .. list-table:: :header-rows: 1 * - - ``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. .. list-table:: :header-rows: 1 * - - ``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``. c Pointer to input/output matrix ``C``. Size of the array must be at least ``ldc`` * ``n``. See :ref:`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 ``alpha`` * op(``A``) * op(``A``)\ :sup:`H` + ``beta`` * ``C``. The imaginary parts of the diagonal elements are set to zero. Return Values ------------- Output event to wait on to ensure computation is complete.