geqrf_batch#

Computes the QR factorizations of a batch of general matrices.

Description

geqrf_batch supports the following precisions.

T

float

double

std::complex<float>

std::complex<double>

geqrf_batch (Buffer Version)#

Description

The buffer version of geqrf_batch supports only the strided API.

Strided API

Syntax

namespace oneapi::mkl::lapack {
  void geqrf_batch(cl::sycl::queue &queue, std::int64_t m, std::int64_t n, cl::sycl::buffer<T> &a, std::int64_t lda, std::int64_t stride_a, cl::sycl::buffer<T> &tau, std::int64_t stride_tau, std::int64_t batch_size, cl::sycl::buffer<T> &scratchpad, std::int64_t scratchpad_size)
}

Input Parameters

queue

Device queue where calculations will be performed.

m

Number of rows in matrices Ai (0m).

n

Number of columns in matrices Ai (0n).

a

Array holding input matrices Ai.

lda

Leading dimension of matrices Ai.

stride_a

Stride between the beginnings of matrices Ai inside the batch array a.

stride_tau

Stride between the beginnings of arrays τi inside the array tau.

batch_size

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 the number of floating point elements of type T. Size should not be less than the value returned by the Strided API of the geqrf_batch_scratchpad_size function.

Output Parameters

a

Factorization data as follows: The elements on and above the diagonal of Ai contain the min(m,n)×n upper trapezoidal matrices Ri (Ri is upper triangular if mn); the elements below the diagonal, with the array τi, contain the orthogonal matrix Qi as a product of min(m,n) elementary reflectors.

tau

Array to store batch of τi, each of size min(m,n), containing scalars that define elementary reflectors for the matrices Qi in its decomposition in a product of elementary reflectors.

Throws

This routine shall throw the following exceptions if the associated condition is detected. An implementation may throw additional implementation-specific exception(s) in case of error conditions not covered here.

oneapi::mkl::lapack::batch_error

oneapi::mkl::unimplemented

oneapi::mkl::unsupported_device

oneapi::mkl::lapack::invalid_argument

The info code of the problem can be obtained by info() method of exception object:

If info = -n, the n-th parameter had an illegal value.

If info equals to value passed as scratchpad size, and detail() returns non zero, then passed scratchpad is of insufficient size, and required size should be not less then value returned by detail() method of exception object.

If info is not zero and detail() returns zero, then there were some errors for some of the problems in the supplied batch and info code contains the number of failed calculations in a batch.

geqrf_batch (USM Version)#

Description

The USM version of geqrf_batch supports the group API and strided API.

Group API

The routine forms the QiRi factorizations of a general m×n matrices Ai, i{1...batch_size}, where batch_size is the sum of all parameter group sizes as provided with group_sizes array. No pivoting is performed during factorization. The routine does not form the matrices Qi explicitly. Instead, Qi is represented as a product of min(m,n) elementary reflectors. Routines are provided to work with Qi in this representation. The total number of problems to solve, batch_size, is a sum of sizes of all of the groups of parameters as provided by group_sizes array.

Syntax

namespace oneapi::mkl::lapack {
  cl::sycl::event geqrf_batch(cl::sycl::queue &queue, std::int64_t *m, std::int64_t *n, T **a, std::int64_t *lda, T **tau, std::int64_t group_count, std::int64_t *group_sizes, T *scratchpad, std::int64_t scratchpad_size, const std::vector<cl::sycl::event> &events = {})
}

Input Parameters

queue

Device queue where calculations will be performed.

m

Array of group_count mg parameters. Each mg specifies the number of rows in matrices Ai from array a, belonging to group g.

n

Array of group_count ng parameters. Each ng specifies the number of columns in matrices Ai from array a, belonging to group g.

a

Array of batch_size pointers to input matrices Ai, each of size ldagng (g is an index of group to which Ai belongs)

lda

Array of group_count ldag parameters, each representing the leading dimensions of input matrices Ai from array a, belonging to group g.

group_count

Specifies the number of groups of parameters. Must be at least 0.

group_sizes

Array of group_count integers. Array element with index g specifies the number of problems to solve for each of the groups of parameters g. So the total number of problems to solve, batch_size, is a sum of all parameter group sizes.

scratchpad

Scratchpad memory to be used by routine for storing intermediate results.

scratchpad_size

Size of scratchpad memory as the number of floating point elements of type T. Size should not be less than the value returned by the Group API of the geqrf_batch_scratchpad_size function.

events

List of events to wait for before starting computation. Defaults to empty list.

Output Parameters

a

Factorization data as follows: The elements on and above the diagonal of Ai contain the min(mg,ng)×ng upper trapezoidal matrices Ri (Ri is upper triangular if mgng); the elements below the diagonal, with the array τi, contain the orthogonal matrix Qi as a product of min(mg,ng) elementary reflectors. Here g is the index of the parameters group corresponding to the i-th decomposition.

tau

Array of pointers to store arrays τi, each of size min(mg,ng), containing scalars that define elementary reflectors for the matrices Qi in its decomposition in a product of elementary reflectors. Here g is the index of the parameters group corresponding to the i-th decomposition.

Return Values

Output event to wait on to ensure computation is complete.

Throws

This routine shall throw the following exceptions if the associated condition is detected. An implementation may throw additional implementation-specific exception(s) in case of error conditions not covered here.

oneapi::mkl::lapack::batch_error

oneapi::mkl::unimplemented

oneapi::mkl::unsupported_device

oneapi::mkl::lapack::invalid_argument

The info code of the problem can be obtained by info() method of exception object:

If info = -n, the n-th parameter had an illegal value.

If info equals to value passed as scratchpad size, and detail() returns non zero, then passed scratchpad is of insufficient size, and required size should be not less then value returned by detail() method of exception object.

If info is not zero and detail() returns zero, then there were some errors for some of the problems in the supplied batch and info code contains the number of failed calculations in a batch.

Strided API

The routine forms the QiRi factorizations of general m×n matrices Ai. No pivoting is performed. The routine does not form the matrices Qi explicitly. Instead, Qi is represented as a product of min(m,n) elementary reflectors. Routines are provided to work with Qi in this representation.

Syntax

namespace oneapi::mkl::lapack {
  sycl::event geqrf_batch(cl::sycl::queue &queue, std::int64_t m, std::int64_t n, T *a, std::int64_t lda, std::int64_t stride_a, T *tau, std::int64_t stride_tau, std::int64_t batch_size, T *scratchpad, std::int64_t scratchpad_size, const std::vector<cl::sycl::event> &events = {})
}

Input Parameters

queue

Device queue where calculations will be performed.

m

Number of rows in matrices Ai (0m).

n

Number of columns in matrices Ai (0n).

a

Array holding input matrices Ai.

lda

Leading dimensions of Ai.

stride_a

Stride between the beginnings of matrices Ai inside the batch array a.

stride_tau

Stride between the beginnings of arrays τi inside the array tau.

batch_size

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 the number of floating point elements of type T. Size should not be less than the value returned by the Strided API of the geqrf_batch_scratchpad_size function.

events

List of events to wait for before starting computation. Defaults to empty list.

Output Parameters

a

Factorization data as follows: The elements on and above the diagonal of Ai contain the min(m,n)×n upper trapezoidal matrices Ri (Ri is upper triangular if mn); the elements below the diagonal, with the array τi, contain the orthogonal matrix Qi as a product of min(m,n) elementary reflectors.

tau

Array to store batch of τi, each of size min(m,n), containing scalars that define elementary reflectors for the matrices Qi in its decomposition in a product of elementary reflectors.

Return Values

Output event to wait on to ensure computation is complete.

Throws

This routine shall throw the following exceptions if the associated condition is detected. An implementation may throw additional implementation-specific exception(s) in case of error conditions not covered here.

oneapi::mkl::lapack::batch_error

oneapi::mkl::unimplemented

oneapi::mkl::unsupported_device

oneapi::mkl::lapack::invalid_argument

The info code of the problem can be obtained by info() method of exception object:

If info = -n, the n-th parameter had an illegal value.

If info equals to value passed as scratchpad size, and detail() returns non zero, then passed scratchpad is of insufficient size, and required size should be not less then value returned by detail() method of exception object.

If info is not zero and detail() returns zero, then there were some errors for some of the problems in the supplied batch and info code contains the number of failed calculations in a batch.

Parent topic: LAPACK-like Extensions Routines