QR factorization

Implementation

Subroutines

Test driver

Example driver

Given \(A \in \mathbb{F}^{m \times n}\), a QR factorization typically computes an implicit unitary matrix \(\hat Q \in \mathbb{F}^{m \times m}\) such that \(\hat R \equiv \hat Q^H A\) is upper trapezoidal. One can then form the thin factors \(Q \in \mathbb{F}^{m \times \mbox{min}(m,n)}\) and \(R \in \mathbb{F}^{\mbox{min}(m,n) \times n}\) by setting \(Q\) and \(R\) to first \(\mbox{min}(m,n)\) columns and rows of \(\hat Q\) and \(\hat R\), respectively. Upon completion \(R\) is stored in the upper trapezoid of \(A\) and the Householder reflectors representing \(\hat Q\) are stored within the columns of the strictly lower trapezoid (this unitary matrix is scaled from the right by a unitary diagonal matrix with entries given by d so that \(R\) has a positive diagonal).

Factorization

C++ API

void QR(Matrix<F> &A, Matrix<F> &t, Matrix<Base<F>> &d)
void QR(AbstractDistMatrix<F> &A, AbstractDistMatrix<F> &t, AbstractDistMatrix<Base<F>> &d)

Overwrite the matrix \(A\) with both \(R\) and the Householder reflectors (and subsequent unitary diagonal matrix defined by the vector, d) representing \(\hat Q\). The scalings for the Householder reflectors are stored in the vector t.

void QR(Matrix<F> &A, Matrix<F> &t, Matrix<Base<F>> &d, Matrix<int> &p, const QRCtrl<Base<F>> &ctrl = QRCtrl<Base<F>>())
void QR(AbstractDistMatrix<F> &A, AbstractDistMatrix<F> &t, AbstractDistMatrix<Base<F>> &d, AbstractDistMatrix<int> &p, const QRCtrl<Base<F>> &ctrl = QRCtrl<Base<F>>())

Overwrite \(A\) with both the \(R\) and (scaled) Householder reflectors from a column-pivoted QR factorization and additionally return the permutation vector, p.

void qr::ExplicitTriang(Matrix<F> &A, const QRCtrl<Base<F>> &ctrl = QRCtrl<Base<F>>())
void qr::ExplicitTriang(AbstractDistMatrix<F> &A, const QRCtrl<Base<F>> &ctrl = QRCtrl<Base<F>>())

Overwrite \(A\) with \(R\).

void qr::ExplicitUnitary(Matrix<F> &A, const QRCtrl<Base<F>> &ctrl = QRCtrl<Base<F>>())
void qr::ExplicitUnitary(AbstractDistMatrix<F> &A, const QRCtrl<Base<F>> &ctrl = QRCtrl<Base<F>>())

Overwrite \(A\) with the orthogonal matrix from its QR factorization (with or without column pivoting).

void qr::Explicit(Matrix<F> &A, Matrix<F> &R, const QRCtrl<Base<F>> &ctrl = QRCtrl<Base<F>>())
void qr::Explicit(AbstractDistMatrix<F> &A, AbstractDistMatrix<F> &R, const QRCtrl<Base<F>> &ctrl = QRCtrl<Base<F>>())

Explicitly return both \(Q\) and \(R\) from the QR factorization.

void qr::Explicit(Matrix<F> &A, Matrix<F> &R, Matrix<Int> &P, const QRCtrl<Base<F>> &ctrl = QRCtrl<Base<F>>())
void qr::Explicit(AbstractDistMatrix<F> &A, AbstractDistMatrix<F> &R, AbstractDistMatrix<int> &P, const QRCtrl<Base<F>> &ctrl = QRCtrl<Base<F>>())

Return representations of all matrices of the pivoted QR factorization. Note that column pivoting is performed regardless of the value of qrCtrl.colPiv.

void qr::Cholesky(Matrix<F> &A, Matrix<F> &R)
void qr::Cholesky(AbstractDistMatrix<F> &A, AbstractDistMatrix<F> &R)

Attempt to perform a QR factorization of a tall-skinny matrix using Cholesky factorization.

qr::TreeData<F> qr::TS(const AbstractDistMatrix<F> &A)

Forms an implicit tall-skinny QR decomposition.

void qr::ExplicitTS(AbstractDistMatrix<F> &A, AbstractDistMatrix<F> &R)

Forms an explicit QR decomposition using a tall-skinny algorithm: A is overwritten with Q.

DistMatrix<F, STAR, STAR> qr::ts::FormR(const AbstractDistMatrix<F> &A, const qr::TreeData<F> &treeData)

Return the R from the QR decomposition.

void qr::ts::FormQ(AbstractDistMatrix<F> &A, qr::TreeData<F> &treeData)

Overwrite A with the Q from the QR decomposition.

template<>
type QRCtrl<Real>
bool colPiv
bool boundRank
Int maxRank
bool adaptive
Real tol
bool alwaysRecomputeNorms
QRCtrl()

Initializes colPiv=false, boundRank=false, maxRank=0, adaptive=false, tol=0, and alwaysRecomputeNorms=false.

template<>
type TreeData<F>
Matrix<F> QR0

Initial QR factorization

Matrix<F> t0

Phases from initial QR factorization

Matrix<Base<F>> d0

Signature (-1,+1) which scales the Householder matrix from the right.

std::vector<Matrix<F>> QRList

Factorizations within reduction tree

std::vector<Matrix<F>> tList

Phases within reduction tree

std::vector<Matrix<Base<F>>> dList

Signatures within reduction tree

C API

ElError ElQR_s(ElMatrix_s A, ElMatrix_s t, ElMatrix_s d)
ElError ElQR_d(ElMatrix_d A, ElMatrix_d t, ElMatrix_d d)
ElError ElQR_c(ElMatrix_c A, ElMatrix_c t, ElMatrix_s d)
ElError ElQR_z(ElMatrix_z A, ElMatrix_z t, ElMatrix_d d)
ElError ElQRDist_s(ElDistMatrix_s A, ElDistMatrix_s t, ElDistMatrix_s d)
ElError ElQRDist_d(ElDistMatrix_d A, ElDistMatrix_d t, ElDistMatrix_d d)
ElError ElQRDist_c(ElDistMatrix_c A, ElDistMatrix_c t, ElDistMatrix_s d)
ElError ElQRDist_z(ElDistMatrix_z A, ElDistMatrix_z t, ElDistMatrix_d d)

Return the packed QR factorization.

ElError ElQRColPiv_s(ElMatrix_s A, ElMatrix_s t, ElMatrix_s d, ElMatrix_i p)
ElError ElQRColPiv_d(ElMatrix_d A, ElMatrix_d t, ElMatrix_d d, ElMatrix_i p)
ElError ElQRColPiv_c(ElMatrix_c A, ElMatrix_c t, ElMatrix_s d, ElMatrix_i p)
ElError ElQRColPiv_z(ElMatrix_z A, ElMatrix_z t, ElMatrix_d d, ElMatrix_i p)
ElError ElQRColPivDist_s(ElDistMatrix_s A, ElDistMatrix_s t, ElDistMatrix_s d, ElDistMatrix_i p)
ElError ElQRColPivDist_d(ElDistMatrix_d A, ElDistMatrix_d t, ElDistMatrix_d d, ElDistMatrix_i p)
ElError ElQRColPivDist_c(ElDistMatrix_c A, ElDistMatrix_c t, ElDistMatrix_s d, ElDistMatrix_i p)
ElError ElQRColPivDist_z(ElDistMatrix_z A, ElDistMatrix_z t, ElDistMatrix_d d, ElDistMatrix_i p)

Return the packed pivoted QR factorization.

ElError ElQRColPivX_s(ElMatrix_s A, ElMatrix_s t, ElMatrix_s d, ElMatrix_i p, ElQRCtrl_s ctrl)
ElError ElQRColPivX_d(ElMatrix_d A, ElMatrix_d t, ElMatrix_d d, ElMatrix_i p, ElQRCtrl_d ctrl)
ElError ElQRColPivX_c(ElMatrix_c A, ElMatrix_c t, ElMatrix_s d, ElMatrix_i p, ElQRCtrl_s ctrl)
ElError ElQRColPivX_z(ElMatrix_z A, ElMatrix_z t, ElMatrix_d d, ElMatrix_i p, ElQRCtrl_d ctrl)
ElError ElQRColPivXDist_s(ElDistMatrix_s A, ElDistMatrix_s t, ElDistMatrix_s d, ElDistMatrix_i p, ElQRCtrl_s ctrl)
ElError ElQRColPivXDist_d(ElDistMatrix_d A, ElDistMatrix_d t, ElDistMatrix_d d, ElDistMatrix_i p, ElQRCtrl_d ctrl)
ElError ElQRColPivXDist_c(ElDistMatrix_c A, ElDistMatrix_c t, ElDistMatrix_s d, ElDistMatrix_i p, ElQRCtrl_s ctrl)
ElError ElQRColPivXDist_z(ElDistMatrix_z A, ElDistMatrix_z t, ElDistMatrix_d d, ElDistMatrix_i p, ElQRCtrl_d ctrl)

Return the packed QR factorization (expert version).

ElError ElQRExplicit_s(ElMatrix_s A, ElMatrix_s R)
ElError ElQRExplicit_d(ElMatrix_d A, ElMatrix_d R)
ElError ElQRExplicit_c(ElMatrix_c A, ElMatrix_c R)
ElError ElQRExplicit_z(ElMatrix_z A, ElMatrix_z R)
ElError ElQRExplicitDist_s(ElDistMatrix_s A, ElDistMatrix_s R)
ElError ElQRExplicitDist_d(ElDistMatrix_d A, ElDistMatrix_d R)
ElError ElQRExplicitDist_c(ElDistMatrix_c A, ElDistMatrix_c R)
ElError ElQRExplicitDist_z(ElDistMatrix_z A, ElDistMatrix_z R)

Return the explicit QR factorization (replace A with Q and return R).

ElError ElQRExplicitColPiv_s(ElMatrix_s A, ElMatrix_s R, ElMatrix_i P)
ElError ElQRExplicitColPiv_d(ElMatrix_d A, ElMatrix_d R, ElMatrix_i P)
ElError ElQRExplicitColPiv_c(ElMatrix_c A, ElMatrix_c R, ElMatrix_i P)
ElError ElQRExplicitColPiv_z(ElMatrix_z A, ElMatrix_z R, ElMatrix_i P)
ElError ElQRExplicitColPivDist_s(ElDistMatrix_s A, ElDistMatrix_s R, ElDistMatrix_i P)
ElError ElQRExplicitColPivDist_d(ElDistMatrix_d A, ElDistMatrix_d R, ElDistMatrix_i P)
ElError ElQRExplicitColPivDist_c(ElDistMatrix_c A, ElDistMatrix_c R, ElDistMatrix_i P)
ElError ElQRExplicitColPivDist_z(ElDistMatrix_z A, ElDistMatrix_z R, ElDistMatrix_i P)

Return the explicit QR factorization with column pivoting (replace A with Q and return R and P).

ElError ElQRExplicitTriang_s(ElMatrix_s A)
ElError ElQRExplicitTriang_d(ElMatrix_d A)
ElError ElQRExplicitTriang_c(ElMatrix_c A)
ElError ElQRExplicitTriang_z(ElMatrix_z A)
ElError ElQRExplicitTriangDist_s(ElDistMatrix_s A)
ElError ElQRExplicitTriangDist_d(ElDistMatrix_d A)
ElError ElQRExplicitTriangDist_c(ElDistMatrix_c A)
ElError ElQRExplicitTriangDist_z(ElDistMatrix_z A)

Return the triangular factor from QR with no pivoting

Note

An expert wrapper which supports column-pivoting is needed.

ElError ElQRExplicitUnitary_s(ElMatrix_s A)
ElError ElQRExplicitUnitary_d(ElMatrix_d A)
ElError ElQRExplicitUnitary_c(ElMatrix_c A)
ElError ElQRExplicitUnitary_z(ElMatrix_z A)
ElError ElQRExplicitUnitaryDist_s(ElDistMatrix_s A)
ElError ElQRExplicitUnitaryDist_d(ElDistMatrix_d A)
ElError ElQRExplicitUnitaryDist_c(ElDistMatrix_c A)
ElError ElQRExplicitUnitaryDist_z(ElDistMatrix_z A)

Return the unitary factor from QR with no pivoting

Note

An expert wrapper which supports column-pivoting is needed.

ElError ElCholeskyQR_s(ElMatrix_s A, ElMatrix_s R)
ElError ElCholeskyQR_d(ElMatrix_d A, ElMatrix_d R)
ElError ElCholeskyQR_c(ElMatrix_c A, ElMatrix_c R)
ElError ElCholeskyQR_z(ElMatrix_z A, ElMatrix_z R)
ElError ElCholeskyQRDist_s(ElDistMatrix_s A, ElDistMatrix_s R)
ElError ElCholeskyQRDist_d(ElDistMatrix_d A, ElDistMatrix_d R)
ElError ElCholeskyQRDist_c(ElDistMatrix_c A, ElDistMatrix_c R)
ElError ElCholeskyQRDist_z(ElDistMatrix_z A, ElDistMatrix_z R)

Attempt to perform a Cholesky-based QR factorization of a tall-skinny matrix.

Apply the factorization to vectors

Applies the implicitly-defined \(Q\) (or its adjoint) stored within A, t, and d from either the left or the right to \(B\).

C++ API

void qr::ApplyQ(LeftOrRight side, Orientation orientation, const Matrix<F> &A, const Matrix<F> &t, const Matrix<Base<F>> &d, Matrix<F> &B)
void qr::ApplyQ(LeftOrRight side, Orientation orientation, const AbstractDistMatrix<F> &A, const AbstractDistMatrix<F> &t, const AbstractDistMatrix<Base<F>> &d, AbstractDistMatrix<F> &B)

C API

ElError ElApplyQAfterQR_s(ElLeftOrRight side, ElOrientation orientation, ElConstMatrix_s A, ElConstMatrix_s t, ElConstMatrix_s d, ElMatrix_s B)
ElError ElApplyQAfterQR_d(ElLeftOrRight side, ElOrientation orientation, ElConstMatrix_d A, ElConstMatrix_d t, ElConstMatrix_d d, ElMatrix_d B)
ElError ElApplyQAfterQR_c(ElLeftOrRight side, ElOrientation orientation, ElConstMatrix_c A, ElConstMatrix_c t, ElConstMatrix_s d, ElMatrix_c B)
ElError ElApplyQAfterQR_z(ElLeftOrRight side, ElOrientation orientation, ElConstMatrix_z A, ElConstMatrix_z t, ElConstMatrix_d d, ElMatrix_z B)
ElError ElApplyQAfterQRDist_s(ElLeftOrRight side, ElOrientation orientation, ElConstDistMatrix_s A, ElConstDistMatrix_s t, ElConstDistMatrix_s d, ElDistMatrix_s B)
ElError ElApplyQAfterQRDist_d(ElLeftOrRight side, ElOrientation orientation, ElConstDistMatrix_d A, ElConstDistMatrix_d t, ElConstDistMatrix_d d, ElDistMatrix_d B)
ElError ElApplyQAfterQRDist_c(ElLeftOrRight side, ElOrientation orientation, ElConstDistMatrix_c A, ElConstDistMatrix_c t, ElConstDistMatrix_s d, ElDistMatrix_c B)
ElError ElApplyQAfterQRDist_z(ElLeftOrRight side, ElOrientation orientation, ElConstDistMatrix_z A, ElConstDistMatrix_z t, ElConstDistMatrix_d d, ElDistMatrix_z B)

Solve linear systems with the factorization

Solves a set of linear systems using an existing packed QR factorization given by \(A\) and the vectors \(t\) and \(d\). \(B\) is the matrix of input vectors and \(X\) is the matrix of solutions.

C++ API

void qr::SolveAfter(Orientation orientation, const Matrix<F> &A, const Matrix<F> &t, const Matrix<Base<F>> &d, const Matrix<F> &B, Matrix<F> &X)
void qr::SolveAfter(Orientation orientation, const AbstractDistMatrix<F> &A, const AbstractDistMatrix<F> &t, const AbstractDistMatrix<Base<F>> &d, const AbstractDistMatrix<F> &B, AbstractDistMatrix<F> &X)

C API

ElError ElSolveAfterQR_s(ElOrientation orientation, ElConstMatrix_s A, ElConstMatrix_s t, ElConstMatrix_s d, ElConstMatrix_s B, ElMatrix_s X)
ElError ElSolveAfterQR_d(ElOrientation orientation, ElConstMatrix_d A, ElConstMatrix_d t, ElConstMatrix_d d, ElConstMatrix_d B, ElMatrix_d X)
ElError ElSolveAfterQR_c(ElOrientation orientation, ElConstMatrix_c A, ElConstMatrix_c t, ElConstMatrix_s d, ElConstMatrix_c B, ElMatrix_c X)
ElError ElSolveAfterQR_z(ElOrientation orientation, ElConstMatrix_z A, ElConstMatrix_z t, ElConstMatrix_d d, ElConstMatrix_z B, ElMatrix_z X)
ElError ElSolveAfterQRDist_s(ElOrientation orientation, ElConstDistMatrix_s A, ElConstDistMatrix_s t, ElConstDistMatrix_s d, ElConstDistMatrix_s B, ElDistMatrix_s X)
ElError ElSolveAfterQRDist_d(ElOrientation orientation, ElConstDistMatrix_d A, ElConstDistMatrix_d t, ElConstDistMatrix_d d, ElConstDistMatrix_d B, ElDistMatrix_d X)
ElError ElSolveAfterQRDist_c(ElOrientation orientation, ElConstDistMatrix_c A, ElConstDistMatrix_c t, ElConstDistMatrix_s d, ElConstDistMatrix_c B, ElDistMatrix_c X)
ElError ElSolveAfterQRDist_z(ElOrientation orientation, ElConstDistMatrix_z A, ElConstDistMatrix_z t, ElConstDistMatrix_d d, ElConstDistMatrix_z B, ElDistMatrix_z X)