LQ Factorization#

GELQF (GEneral LQ Factorization) function computes batched LQ factorization of a general matrix:

\[A = L * Q,\]

using Householder reflection transformations, where:

  • A is a batched M x N matrix, with leading dimension lda >= M if matrix A is in column-major layout, or lda >= N if matrix A is row-major,

  • Q is an unitary N x N matrix, and

  • L is a lower triangular matrix if M <= N, or lower trapezoidal matrix if M > N.

The LQ factorization of A is essentially the same as the QR factorization of \(A^T\), or \(A^H\) for complex data type.

cuSolverDx gelqf device functions are (see Execution Methods):

__device__ void execute(data_type* A, data_type* tau);
// with runtime lda
__device__ void execute(data_type* A, const unsigned int lda, data_type* tau);

After the function, the lower triangular or lower trapezoidal part of input A, including diagonal elements, is replaced by the matrix L.

Matrix Q is not explicitly formed. The elements above the diagonal of A, with the array tau, represent the matrix Q as a product of min(M, N) Householder vectors:

\[Q = H(min(M, N) -1 )^H * . . . * H(1)^H * H(0)^H.\]

Each Householder vector has the form \(H(i) = I - tau[i] * v * v^H\), where:

  • tau is an array of size min(M, N) for each batch, and

  • v is a vector of size N for each batch, with v[0 : i - 1] = 0, v[i] = 1, and conjugate(v[i+1 : N]) is stored on exit in A[i, i+1 : N].

The functions support A being either column- or row-major memory layout, see Arrangement operator.