gemv_batch
Computes a group of
gemv
operations.Description
The
gemv_batch
routines are batched versions of gemv, performing multiple gemv
operations in a single call.
Each gemv
operations perform a scalar-matrix-vector product and add the result to a scalar-vector product.gemv_batch
supports the following precisions:T |
---|
float |
double |
std::complex<float> |
std::complex<double> |
gemv_batch (Buffer Version)
Buffer version of
gemv_batch
supports only strided API.Strided API
Strided API operation is defined as:
for i = 0 … batch_size – 1
A is a matrix at offset i * stridea in a.
X and Y are vctors at offset i * stridex, i * stridey, in x and y.
Y = alpha * op(A) * X + beta * Y
end for
where:
- op(A) is one of op(A) =A, or op(A) =AT, or op(A) =AH
- alphaandbetaare scalars
- Ais matrix andXandYare vectors
For strided API,
x
and y
buffers contain all the input vectors. The stride between vectors is either given by the stride parameters. Total number of vectors in x
and y
buffers is given by batch_size
parameter.Syntax
namespace oneapi::mkl::blas::column_major {
void gemv_batch(sycl::queue &queue,
oneapi::mkl::transpose trans,
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> &x,
std::int64_t incx,
std::int64_t stridex,
T beta,
sycl::buffer<T,1> &y,
std::int64_t incy,
std::int64_t stridey,
std::int64_t batch_size)
}
namespace oneapi::mkl::blas::row_major {
void gemv_batch(sycl::queue &queue,
oneapi::mkl::transpose trans,
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> &x,
std::int64_t incx,
std::int64_t stridex,
T beta,
sycl::buffer<T,1> &y,
std::int64_t incy,
std::int64_t stridey,
std::int64_t batch_size)
}
Input Parameters
- queue
- The queue where the routine should be executed.
- trans
- Specifies op(A), the transposition operation applied to matricesA. See Data Types for more details.
- m
- Number of rows of matrices op(A). Must be at least zero.
- n
- Number of columns of matrices op(A). Must be at least zero.
- alpha
- Scaling factor for matrix-vector product.
- a
- The buffer holding input matricesA. Size of the buffer must be at leaststridea*batch_size.
- lda
- Leading dimension of matricesA. Must be positive and at leastmif column major layout or at leastnif row major layout is used.
- stridea
- Stride between two consecutiveAmatrices.
- x
- Buffer holding input vectorsX. Size of the buffer must be at leaststridex*batch_size.
- incx
- Stride between two consecutive elements ofXvectors.
- stridex
- Stride between two consecutiveXvectors. Must be at least zero.
- beta
- Scaling factor for vectorsY.
- y
- Buffer holding input/output vectorsY. Size of the buffer must be at leaststridey*batch_size.
- incy
- Stride between two consecutive elements ofYvectors.
- stridey
- Stride between two consecutiveYvectors. Must be at least zero.
- batch_size
- Number ofgemvcomputations to perform. Must be at least zero.
Output Parameters
- y
- Output buffer overwritten bybatch_sizegemvoperations of the formalpha* op(A) *X+beta*Y.
gemv_batch (USM Version)
USM version of
gemv_batch
supports group API 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 is an m x n matrix in a[idx]
X and Y are vectors in x[idx] and y[idx]
Y = alpha[i] * op(A) * X + beta[i] * Y
idx = idx + 1
end for
end for
where:
- op(A) is one of op(A) =A, or op(A) =AT, or op(A) =AH
- alphaandbetaare scalars
- Ais matrix andXandYare vectors
For group API,
x
and y
arrays contain the pointers for all the input vectors.
a
array contains the pointers to all input matrices.
The total number of vectors in x
and y
and matrices in a
is given by:
Syntax
namespace oneapi::mkl::blas::column_major {
sycl::event gemv_batch(sycl::queue &queue,
oneapi::mkl::transpose *trans,
std::int64_t *m,
std::int64_t *n,
T *alpha,
const T **a,
std::int64_t *lda,
const T **x,
std::int64_t *incx,
T *beta,
T **y,
std::int64_t *incy,
std::int64_t group_count,
std::int64_t *group_size,
const std::vector<sycl::event> &dependencies = {})
}
namespace oneapi::mkl::blas::row_major {
sycl::event gemv_batch(sycl::queue &queue,
oneapi::mkl::transpose *trans,
std::int64_t *m,
std::int64_t *n,
T *alpha,
const T **a,
std::int64_t *lda,
const T **x,
std::int64_t *incx,
T *beta,
T **y,
std::int64_t *incy,
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.
- trans
- Array ofgroup_countoneapi::mkl::transposevalues.transa[i]specifies op(A), the transposition operation applied to matricesAin groupi. See Data Types for more details.
- m
- Array ofgroup_countintegers.m[i]specifies number of rows of matrices op(A) in groupi. All entries must be at least zero.
- n
- Array ofgroup_countintegers.n[i]specifies number of columns of matrices op(A) in groupi. All entries must be at least zero.
- alpha
- Array ofgroup_countscalar elements.alpha[i]specifies scaling factor for matrix-vector products in groupi.
- a
- lda
- Array ofgroup_countintegers.lda[i]specifies leading dimension of matricesAin groupi. Must be positive and at leastm[i]if column major layout or at leastn[i]if row major layout is used.
- x
- incx
- Array ofgroup_countintegers.incx[i]specifies stride of vectorsXin groupi.
- beta
- Array ofgroup_countscalar elements.beta[i]specifies scaling factor for vectorsYin groupi.
- y
- Array oftotal_batch_countpointers for input/output vectorsY. See Matrix Storage for more details.
- incy
- Array ofgroup_countintegers.incy[i]specifies stride of vectorsYin groupi.
- group_count
- Number of groups. Must be at least zero.
- group_size
- Array ofgroup_countintegers.group_size[i]specifies the number ofgemvoperations 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
- y
- Array of pointers to output vectorsYoverwritten bytotal_batch_countgemvoperations of the formalpha* op(A) *X+beta*Y.
Return Values
Output event to wait on to ensure computation is complete.
Examples
An example of how to use USM version of
gemv_batch
can be found in oneMKL installation directory, under:examples/dpcpp/blas/source/gemv_batch_usm.cpp
Strided API
Strided API operation is defined as:
for i = 0 … batch_size – 1
A is a matrix at offset i * stridea in a.
X and Y are vctors at offset i * stridex, i * stridey, in x and y.
Y = alpha * op(A) * X + beta * Y
end for
where:
- op(A) is one of op(A) =A, or op(A) =AT, or op(A) =AH
- alphaandbetaare scalars
- Ais matrix andXandYare vectors
For strided API,
x
and y
arrays contain all the input vectors. The stride between vectors is given by the stride parameters. Total number of vectors in x
and y
arrays is given by batch_size
parameter.Syntax
namespace oneapi::mkl::blas::column_major {
sycl::event gemv_batch(sycl::queue &queue,
oneapi::mkl::transpose trans,
std::int64_t m,
std::int64_t n,
T alpha,
const T *a,
std::int64_t lda,
std::int64_t stridea,
const T *x,
std::int64_t incx,
std::int64_t stridex,
T beta,
T *y,
std::int64_t incy,
std::int64_t stridey,
std::int64_t batch_size)
}
namespace oneapi::mkl::blas::row_major {
sycl::event gemv_batch(sycl::queue &queue,
oneapi::mkl::transpose trans,
std::int64_t m,
std::int64_t n,
T alpha,
const T *a,
std::int64_t lda,
std::int64_t stridea,
const T *x,
std::int64_t incx,
std::int64_t stridex,
T beta,
T *y,
std::int64_t incy,
std::int64_t stridey,
std::int64_t batch_size)
}
Input Parameters
- queue
- The queue where the routine should be executed.
- trans
- Specifies op(A), the transposition operation applied to matricesA. See Data Types for more details.
- m
- Number of rows of matrices op(A). Must be at least zero.
- n
- Number of columns of matrices op(A). Must be at least zero.
- alpha
- Scaling factor for matrix-vector product.
- a
- Pointer to input matricesA. Size of the array must be at leaststridea*batch_size.
- lda
- Leading dimension of matricesA. Must be positive and at leastmif column major layout or at leastnif row major layout is used.
- stridea
- Stride between two consecutiveAmatrices.
- x
- Pointer to input vectorsX. Size of the array must be at leaststridex*batch_size.
- incx
- Stride between two consecutive elements ofXvectors.
- stridex
- Stride between two consecutiveXvectors. Must be at least zero.
- beta
- Scaling factor for vectorsY.
- y
- Pointer to input/output vectorsY. Size of the array must be at leaststridey*batch_size.
- incy
- Stride between two consecutive elements ofYvectors.
- stridey
- Stride between two consecutiveYvectors. Must be at least zero.
- batch_size
- Number ofgemvcomputations to perform. Must be at least zero.
- dependencies
- List of events to wait for before starting computation, if any. If omitted, defaults to no dependencies.
Output Parameters
- y
- Pointer to output vectorsYoverwritten bybatch_sizegemvoperations of the formalpha* op(A) *X+beta*Y.
Return Values
Output event to wait on to ensure computation is complete.