gesvda_batch (USM Strided Version)#
Computes the truncated SVD factorizations for a batch of general matrices. This routine belongs to the oneapi::mkl::lapack
namespace.
Description#
The routine computes the truncated SVD decomposition (or a reduced rank
approximation) for a batch of general matrices A
i, as:
A
i = U
i * S
i * V
iT
where U
i and V
i are orthogonal matrices, S
i
is a diagonal matrix with singular values on the diagonal. Singular
values are nonnegative and listed in decreasing order. A truncated SVD
of a given matrix will produce matrices with the specified number of
columns, where the number of columns is defined by the user or
determined at runtime with the help of the user-defined tolerance
threshold.
An approximation of each matrix can be also obtained as a product of two low-rank matrices (e.g. low rank product):
A
i = P
i * Q
i
where
P
i = U
i * S
i,
Q
i = V
iT
if the number of rows is greater or equal to the number of columns,
P
i = U
i,
Q
i = S
i` *V
iT
otherwise.
The routine provides 3 possible ways for computing truncated SVD:
Compute truncated SVD with the help of the input array
rank
whererank[i]
specifies the number of singular values and vectors to be computed in parameterU
i,V
i andS
i for each matrixA
i.Compute truncated SVD using a tolerance threshold. While computing SVD, singular values that are less than the user-defined tolerance are treated as zero, and they are not computed but set to zero.
Compute truncated SVD using the effective rank. The effective rank is determined by treating as zero those singular values that are less than the user-defined tolerance threshold times the largest singular value.
The routine can be also used for computing singular values only.
API#
Syntax#
namespace oneapi::mkl::lapack {
void gesvda_batch(sycl::queue &queue,
int64_t *iparm,
int64_t *irank,
int64_t m,
int64_t n,
T *a,
int64_t lda,
int64_t stride_a,
RealT *s,
int64_t stride_s,
T * u,
int64_t ldu,
int64_t stride_u,
T * vt,
int64_t ldvt,
int64_t stride_vt,
RealT tolerance,
RealT *residual,
int64_t batch_size,
T *scratchpad,
int64_t scratchpad_size,
const std::vector<sycl::event> &events = {})
}
gesvda_batch
supports the following precisions and devices:
T |
Devices supported |
---|---|
|
CPU, GPU* |
|
CPU, GPU* |
|
CPU, GPU* |
|
CPU, GPU* |
*Interface support only; all computations are performed on the CPU.
Input Parameters#
- queue
Device queue where calculations will be performed.
- iparm
Array of dimension 16 specifying options to compute truncated SVD. It also specifies the type of returned SVD decomposition form. The table below describes all individual components of the
iparm
parameter. Default values are denoted with an asterisk (*). NOTE:iparm[4]-iparm[15]
are reserved for future use.
Component description |
Values |
---|---|
|
If If If If |
|
If If |
|
If
is computed. If
is computed. |
|
If If |
- irank
If
iparm[0]=0
oriparm[0]=-1
, elementirank[i]
specifies the number of singular values and/or singular vectors to be computed inU
i,V
iTandS
i for each matrixA
i.- m
The number of rows in the matrices
A
i.- n
The number of columns in the matrices
A
i.- a
Array holding input matrices
A
i.- lda
The leading dimension of a. Must be at least
max(1, m)
.- stride_a
The stride between the beginnings of matrices
A
i inside the batch array a. Must be at leastmax(1, lda * n)
.- stride_s
The stride between the beginnings of arrays
S
i inside the array s. Must be at leastmin(m,n)
.- ldu
The leading dimension of
U
i. Must be at leastmax(1, m)
.- stride_u
The stride between the beginnings of matrices
U
i inside the batch array u. Must be at leastmax(1, ldu * m)
.- ldvt
The leading dimension of
V
iT . Must be at leastmax(1, n)
.- stride_vt
The stride between the beginnings of matrices
V
iTinside the batch array vt. Must be at leastmax(1, ldvt * n)
.- tolerance
Specifies the tolerance threshold. It is only used for computing truncated SVD in the cases of
iparm[0]=1
andiparm[0]=2
. Not used otherwise.- batch_size
The number of problems in a batch.
- scratchpad
Scratchpad memory to be used by routine for storing intermediate results.
- scratchpad_size
Size of scratchpad memory as a number of floating point elements of type T. Size should not be less than the value returned by gesvda_batch_scratchpad_size (Strided Version) function.
Output Parameters#
- irank
If
iparm[0]=-1
oriparm[0]=0
, elementirank[i]
is the number of computed singular values and/or singular vectors for matrixA
i.- a
Unchanged on exit, if the residual vector is not required. Otherwise contains the residual matrix
A
i -U
i *S
i *V
iTifiparm[2]=0
,A
i -P
i *Q
i ifiparm[2]=1
.- s
Array of size
min(m,n)*batch_size
to store batch of singular valuesS
i.- u
u is array of size at least
stride_u*batch_size
to store batch ofU
i ifiparm[2]=0
, or to store batch ofP
i ifiparm[2]=1
.- vt
vt is array of size at least
stride_vt*batch_size
to store batch ofV
iTifiparm[2]=0
, or to store batch ofQ
i ifiparm[2]=1
.- residual
Array of dimension
batch_size
. Ifiparm[3]=1
,residual[i]
is the Frobenius norm of the matrix:||
A
i -U
i *S
i *V
iT ||if
iparm[2]=0
, and||
A
i -P
i *Q
i ||if
iparm[2]=1
.
Exceptions#
Exception |
Description |
---|---|
|
This exception is thrown when problems occur during calculations. You can obtain the info code of the problem using the info() method of the exception object: If If If If If |
Return Values#
Output event to wait on to ensure computation is complete.