trsm_batch
Computes a group of
trsm
operations.Description
The
trsm_batch
routines are batched versions of trsm, performing multiple trsm
operations in a single call. Each trsm
solves an equation of the form op(A) * X = alpha * B or X * op(A) = alpha * B.trsm_batch
supports the following precisions:T |
---|
float |
double |
std::complex<float> |
std::complex<double> |
trsm_batch (Buffer Version)
Buffer version of
trsm_batch
supports only strided API.Strided API
Strided API operation is defined as:
for i = 0 … batch_size – 1
A and B are matrices at offset i * stridea and i * strideb in a and b.
if (left_right == side::left) then
compute X such that op(A) * X = alpha * B
else
compute X such that X * op(A) = alpha * B
B = X
end for
where:
- op(A) is one of op(A) =A, or op(A) =AT, or op(A) =AH
- alphais a scalar
- Ais eithermxmornxntriangular matrix
- BandXaremxngeneral matrices
On return, matrix
B
is overwritten by solution matrix X
.For strided API,
a
and b
buffers contains all the input matrices. The stride between matrices is given by the stride parameters. Total number of matrices in a
and b
buffers is given by batch_size
parameter.Syntax
namespace oneapi::mkl::blas::column_major {
void trsm_batch(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<T,1> &a,
std::int64_t lda,
std::int64_t stridea,
sycl::buffer<T,1> &b,
std::int64_t ldb,
std::int64_t strideb,
std::int64_t batch_size)
}
namespace oneapi::mkl::blas::row_major {
void trsm_batch(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<T,1> &a,
std::int64_t lda,
std::int64_t stridea,
sycl::buffer<T,1> &b,
std::int64_t ldb,
std::int64_t strideb,
std::int64_t batch_size)
}
Input Parameters
- queue
- The queue where the routine should be executed.
- left_right
- Specifies whether matricesAare on the left side or right side of the multiplication. See Data Types for more details.
- upper_lower
- Specifies whether matricesAare upper or lower triangular. See Data Types for more details.
- trans
- unit_diag
- Specifies whether matricesAare unit triangular or not. See Data Types for more details.
- m
- Number of rows of matricesB. Must be at least zero.
- n
- Number of columns of matricesB. Must be at least zero.
- alpha
- Scaling factor for the solution.
- a
- Buffer holding input matriceesA. Size of the buffer must be at leaststridea*batch_size.
- lda
- Leading dimension of matricesA. Must be at leastmifleft_right=side::leftor at leastnifleft_right=side::right. Must be positive.
- stridea
- Stride between two consecutiveAmatrices.
- b
- Buffer holding input/output matricesB. Size of the buffer must be at leaststrideb*batch_size.
- ldb
- Leading dimension of matricesB. Must be at leastmif column major layout or at leastnif row major layout is used. Must be positive.
- strideb
- Stride between two consecutiveBmatrices.
- batch_size
- Specifies number of triangular linear systems to solve.
Output Parameters
- b
- Output buffer overwritten bybatch_sizesolution matricesX.
If
alpha
= 0, matrices B
are set to zero, and A
and B
do not need to be initialized before calling trsm_batch
.. trsm_batch (USM Version)
USM version of
trsm_batch
supports group API and strided API.Group API
Group API operation is defined as:
idx = 0
for i = 0 … group_count – 1
for j = 0 … group_size – 1
A and B are matrices in a[idx] and b[idx]
if (left_right == side::left) then
compute X such that op(A) * X = alpha[i] * B
else
compute X such that X * op(A) = alpha[i] * B
end if
B = X
idx = idx + 1
end for
end for
where:
- op(A) is one of op(A) =A, or op(A) =AT, or op(A) =AH
- alphais a scalar
- Ais eithermxmornxntriangular matrix
- BandXaremxngeneral matrices
On return, matrix
B
is overwritten by solution matrix X
.For group API,
a
and b
arrays contain the pointers for all the input matrices.
The total number of matrices in a
and b
are given by:
Syntax
namespace oneapi::mkl::blas::column_major {
sycl::event trsm_batch(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,
std::int64_t group_count,
std::int64_t *group_size,
const std::vector<sycl::event> &dependencies = {})
}
namespace oneapi::mkl::blas::row_major {
sycl::event trsm_batch(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,
std::int64_t group_count,
std::int64_t *group_size,
const std::vector<sycl::event> &dependencies = {})
}
Input Parameters
- queue
- The queue where the routine should be executed.
- left_right
- Array ofgroup_countoneapi::mkl::sidevalues.left_right[i]specifies whether matricesAare on the left side or right side of the multiplication in groupi. See Data Types for more details.
- upper_lower
- Array ofgroup_countoneapi::mkl::uplovalues.upper_lower[i]specifies whether matricesAare upper or lower triangular in groupi. See Data Types for more details.
- trans
- Array ofgroup_countoneapi::mkl::transposevalues.trans[i]specifies op(A), transposition operation applied to matricesAin each groupi. See Data Types for more details.
- unit_diag
- Array ofgroup_countoneapi::mkl::diagvalues.unit_diag[i]specifies whether matricesAare unit triangular or not. See Data Types for more details.
- m
- Array ofgroup_countintegers.m[i]specifies number of rows of matricesBin groupi. All entries must be at least zero.
- n
- Array ofgroup_countintegers.n[i]specifies number of columns of matricesBin groupi. All entries must be at least zero.
- alpha
- Array ofgroup_countscalar elements.alpha[i]specifies scaling factors for the solutions in groupi.
- a
- lda
- Array ofgroup_countintegers.lda[i]specifies leading dimension of matricesAin groupi. Must be at leastm[i]ifleft_right[i]=side::leftor at leastn[i]ifleft_right[i]=side::right. All entries must be positive.
- b
- Array oftotal_batch_countpointers for input/output matricesB. See Matrix Storage for more details.
- ldb
- Array ofgroup_countintegers.ldb[i]specifies leading dimension of matricesBin groupi. Must be at leastm[i]if column major layout or at leastn[i]if row major layout is used. All entries must be positive.
- group_count
- Number of groups. Must be at least zero.
- group_size
- Array ofgroup_countintegers.group_size[i]specifies the number oftrsmoperations in groupi. Each element ingroup_sizemust be at least zero.
- dependencies
- List of events to wait for before starting computation, if any. If omitted, defaults to no dependencies.
Output Parameters
- b
- Array of pointers to output matricesBoverwritten bytotal_batch_countsolution matricesX.
If
alpha
= 0, matrices B
are set to zero, and A
and B
do not need to be initialized before calling trsm_batch
.. Return Values
Output event to wait on to ensure computation is complete.
Strided API
Strided API operation is defined as:
for i = 0 … batch_size – 1
A and B are matrices at offset i * stridea and i * strideb in a and b.
if (left_right == side::left) then
compute X such that op(A) * X = alpha * B
else
compute X such that X * op(A) = alpha * B
B = X
end for
where:
- op(A) is one of op(A) =A, or op(A) =AT, or op(A) =AH
- alphais a scalar
- Ais eithermxmornxntriangular matrix
- BandXaremxngeneral matrices
On return, matrix
B
is overwritten by solution matrix X
.For strided API,
a
and b
arrays contain all the input matrices. The stride between matrices is given by the stride parameters. Total number of matrices in a
and b
arrays is given by batch_size
parameter.Syntax
namespace oneapi::mkl::blas::column_major {
sycl::event trsm_batch(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,
std::int64_t stridea,
T *b,
std::int64_t ldb,
std::int64_t strideb,
std::int64_t batch_size,
const std::vector<sycl::event> &dependencies = {})
}
namespace oneapi::mkl::blas::row_major {
sycl::event trsm_batch(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,
std::int64_t stridea,
T *b,
std::int64_t ldb,
std::int64_t strideb,
std::int64_t batch_size,
const std::vector<sycl::event> &dependencies = {})
}
Input Parameters
- queue
- The queue where the routine should be executed.
- left_right
- Specifies whether matricesAare on the left side or right side of the multiplication. See Data Types for more details.
- upper_lower
- Specifies whether matricesAare upper or lower triangular. See Data Types for more details.
- trans
- unit_diag
- Specifies whether matricesAare unit triangular or not. See Data Types for more details.
- m
- Number of rows of matricesB. Must be at least zero.
- n
- Number of columns of matricesB. Must be at least zero.
- alpha
- Scaling factor for the solution.
- a
- Pointer to input matriceesA. Size of the array must be at leaststridea*batch_size.
- lda
- Leading dimension of matricesA. Must be at leastmifleft_right=side::leftor at leastnifleft_right=side::right. Must be positive.
- stridea
- Stride between two consecutiveAmatrices.
- b
- Pointer to input/output matricesB. Size of the array must be at leaststrideb*batch_size.
- ldb
- Leading dimension of matricesB. Must be at leastmif column major layout or at leastnif row major layout is used. Must be positive.
- strideb
- Stride between two consecutiveBmatrices.
- batch_size
- Specifies number of triangular linear systems to solve.
- dependencies
- List of events to wait for before starting computation, if any. If omitted, defaults to no dependencies.
Output Parameters
- b
- Pointer to output matrixBoverwritten bybatch_sizesolution matricesX.
If
alpha
= 0, matrices B
are set to zero, and A
and B
do not need to be initialized before calling trsm_batch
.. Return Values
Output event to wait on to ensure computation is complete.