How to find (Q, R ) from SuiteSparseQR_factorization object?

Question

In C++ interface of SuiteSparse, I can use

SuiteSparseQR_factorization <double> *QR;
QR = SuiteSparseQR_factorize(A) ;

to calculate QR decomposition of matrix A so that I can reuse QR for further calculation. But I wonder can I get the real Q,R directly from this QR object?

score 0 · Answer 1 · answered Jun 07 '20 at 10:14

SuiteSparse is awesome, but the interface can be confusing. Unfortunately, the methods that involve the SuiteSparseQR_factorization struct, which appear to be the most convenient, haven't worked so well for me in practice. For instance, using SuiteSparseQR_factorize and then SuiteSparseQR_qmult with a sparse matrix input argument actually converts it to a dense matrix first, which seems completely unnecessary!

Instead, use

template <typename Entry> SuiteSparse_long SuiteSparseQR
(
    // inputs, not modified
    int ordering,           // all, except 3:given treated as 0:fixed
    double tol,             // only accept singletons above tol

    SuiteSparse_long econ,  // number of rows of C and R to return; a value
                            // less than the rank r of A is treated as r, and
                            // a value greater than m is treated as m.

    int getCTX,             // if 0: return Z = C of size econ-by-bncols
                            // if 1: return Z = C' of size bncols-by-econ
                            // if 2: return Z = X of size econ-by-bncols

    cholmod_sparse *A,      // m-by-n sparse matrix

    // B is either sparse or dense.  If Bsparse is non-NULL, B is sparse and
    // Bdense is ignored.  If Bsparse is NULL and Bdense is non-NULL, then B is
    // dense.  B is not present if both are NULL.
    cholmod_sparse *Bsparse,
    cholmod_dense *Bdense,

    // output arrays, neither allocated nor defined on input.

    // Z is the matrix C, C', or X
    cholmod_sparse **Zsparse,
    cholmod_dense  **Zdense,
    cholmod_sparse **R,     // the R factor
    SuiteSparse_long **E,   // size n; fill-reducing ordering of A.
    cholmod_sparse **H,     // the Householder vectors (m-by-nh)
    SuiteSparse_long **HPinv,// size m; row permutation for H
    cholmod_dense **HTau,   // size nh, Householder coefficients

    // workspace and parameters
    cholmod_common *cc
) ;

This method will perform the factorization and then, optionally, output (among other things) R, the matrix product Z = Q^T * B (or its transpose -- B^T * Q), or the solution of a linear system. To get Q, define B as the identity matrix. Here's an example to get Q and R.

cholmod_common Common, * cc;
cc = &Common;
cholmod_l_start(cc);
cholmod_sparse *A;//assume you have already defined this
int ordering = SPQR_ORDERING_BEST;
double tol = 0;
Long econ = A->nrow;
int getCTX = 1;// Z = (Q^T * B)^T = B^T * Q
cholmod_sparse *B = cholmod_l_speye(A->nrow, A->nrow, CHOLMOD_REAL, cc);//the identity matrix
cholmod_sparse *Q, *R;//output pointers to the Q and R sparse matrices

SuiteSparseQR<double>(ordering, tol, econ, getCTX, A, B, NULL, &Q, NULL, &R, NULL, NULL, NULL, NULL, cc);

If you want any of the other outputs to perform subsequent operations without the use of an explicitly formed Q and/or R, then you need to substitute the NULL's for additional pointers and then make calls to SuiteSparseQR_qmult.

How to find (Q, R ) from SuiteSparseQR_factorization object?

1 Answers1