.. _onemkl_blas_trsm: trsm ==== Solves a triangular matrix equation (forward or backward solve). Description *********** The ``trsm`` routines solve triangular matrix equations where one of the matrices in the multiplication is triangular. The argument ``left_right`` determines if the triangular matrix, ``A``, is on the left of the multiplication (``left_right`` = ``side::left``) or on the right (``left_right`` = ``side::right``). The operation is defined as: If (``left_right`` = ``side::left``), .. math:: op(A)*X = alpha*B If (``left_right`` = ``side::right``), .. math:: X*op(A) = alpha*B where: - op(``A``) is one of op(``A``) = ``A``, or op(``A``) = ``A``\ :sup:`T`, or op(``A``) = ``A``\ :sup:`H` - ``alpha`` is a scalar - ``A`` is either ``m`` x ``m`` or ``n`` x ``n`` triangular matrix - ``B`` and ``X`` are ``m`` x ``n`` general matrices On return, matrix ``B`` is overwritten by solution matrix ``X``. ``trsm`` supports the following precisions: .. list-table:: :header-rows: 1 * - T * - ``float`` * - ``double`` * - ``std::complex`` * - ``std::complex`` trsm (Buffer Version) ********************* Syntax ------ .. code-block:: cpp namespace oneapi::mkl::blas::column_major { void trsm(sycl::queue &queue, oneapi::mkl::side left_right, oneapi::mkl::uplo upper_lower, oneapi::mkl::transpose transa, oneapi::mkl::diag unit_diag, std::int64_t m, std::int64_t n, T alpha, sycl::buffer &a, std::int64_t lda, sycl::buffer &b, std::int64_t ldb) } .. code-block:: cpp namespace oneapi::mkl::blas::row_major { void trsm(sycl::queue &queue, oneapi::mkl::side left_right, oneapi::mkl::uplo upper_lower, oneapi::mkl::transpose trans, oneapi::mkl::diag unit_diag, std::int64_t m, std::int64_t n, T alpha, sycl::buffer &a, std::int64_t lda, sycl::buffer &b, std::int64_t ldb) } Input Parameters ---------------- queue The queue where the routine should be executed. left_right Specifies whether matrix ``A`` is on the left side or right side of the multiplication. See :ref:`data-types` for more details. upper_lower Specifies whether matrix ``A`` is upper or lower triangular. See :ref:`data-types` for more details. trans Specifies op(``A``), the transposition operation applied to matrix ``A``. See :ref:`data-types` for more details. unit_diag Specifies whether matrix ``A`` is unit triangular or not. See :ref:`data-types` for more details. m Number of rows of matrix ``B``. Must be at least zero. n Number of columns of matrix ``B``. Must be at least zero. alpha Scaling factor for the solution. a Buffer holding input matrix ``A``. Size of the buffer must be at least ``lda`` * ``m`` if ``left_right`` = ``side::left`` or ``lda`` * ``n`` if ``left_right`` = ``side::right``. See :ref:`matrix-storage` for more details. lda Leading dimension of matrix ``A``. Must be at least ``m`` if ``left_right`` = ``side::left`` or at least ``n`` if ``left_right`` = ``side::right``. Must be positive. b Buffer holding input matrix ``B``. Size of the buffer must be at least ``ldb`` * ``n`` if column major layout or at least ``ldb`` * ``m`` if row major layout is used. See :ref:`matrix-storage` for more details. ldb Leading dimension of matrix ``B``. Must be at least ``m`` if column major layout or at least ``n`` if row major layout is used. Must be positive. Output Parameters ----------------- b Output buffer overwritten by solution matrix ``X``. .. note:: If ``alpha`` = 0, matrix ``B`` is set to zero, and ``A`` and ``B`` do not need to be initialized before calling ``trsm``. trsm (Buffer Version) ********************* Syntax ------ .. code-block:: cpp namespace oneapi::mkl::blas::column_major { sycl::event trsm(sycl::queue &queue, oneapi::mkl::side left_right, oneapi::mkl::uplo upper_lower, oneapi::mkl::transpose trans, oneapi::mkl::diag unit_diag, std::int64_t m, std::int64_t n, T alpha, const T* a, std::int64_t lda, T* b, std::int64_t ldb, const std::vector &dependencies = {}) } .. code-block:: cpp namespace oneapi::mkl::blas::row_major { sycl::event trsm(sycl::queue &queue, oneapi::mkl::side left_right, oneapi::mkl::uplo upper_lower, oneapi::mkl::transpose trans, oneapi::mkl::diag unit_diag, std::int64_t m, std::int64_t n, T alpha, const T* a, std::int64_t lda, T* b, std::int64_t ldb, const std::vector &dependencies = {}) } Input Parameters ---------------- queue The queue where the routine should be executed. left_right Specifies whether matrix ``A`` is on the left side or right side of the multiplication. See :ref:`data-types` for more details. upper_lower Specifies whether matrix ``A`` is upper or lower triangular. See :ref:`data-types` for more details. trans Specifies op(``A``), the transposition operation applied to matrix ``A``. See :ref:`data-types` for more details. unit_diag Specifies whether matrix ``A`` is unit triangular or not. See :ref:`data-types` for more details. m Number of rows of matrix ``B``. Must be at least zero. n Number of columns of matrix ``B``. Must be at least zero. alpha Scaling factor for the solution. a Pointer to input matrix ``A``. Size of the array must be at least ``lda`` * ``m`` if ``left_right`` = ``side::left`` or ``lda`` * ``n`` if ``left_right`` = ``side::right``. See :ref:`matrix-storage` for more details. lda Leading dimension of matrix ``A``. Must be at least ``m`` if ``left_right`` = ``side::left`` or at least ``n`` if ``left_right`` = ``side::right``. Must be positive. b Pointer to input matrix ``B``. Size of the array must be at least ``ldb`` * ``n`` if column major layout or at least ``ldb`` * ``m`` if row major layout is used. See :ref:`matrix-storage` for more details. ldb Leading dimension of matrix ``B``. Must be at least ``m`` if column major layout or at least ``n`` if row major layout is used. Must be positive. dependencies List of events to wait for before starting computation, if any. If omitted, defaults to no dependencies. Output Parameters ----------------- b Pointer to output matrix overwritten by solution matrix ``X``. .. note:: If ``alpha`` = 0, matrix ``B`` is set to zero, and ``A`` and ``B`` do not need to be initialized before calling ``trsm``. Return Values ------------- Output event to wait on to ensure computation is complete.