Factorizations¶
Cholesky factorization¶
It is well-known that Hermitian positive-definite (HPD) matrices can be decomposed into the form \(A = L L^H\) or \(A = U^H U\), where \(L=U^H\) is lower triangular, and Cholesky factorization provides such an \(L\) (or \(U\)) given an HPD \(A\). If \(A\) is found to be numerically indefinite, then a NonHPDMatrixException will be thrown.
- void Cholesky(UpperOrLower uplo, Matrix<F>& A)¶
- void Cholesky(UpperOrLower uplo, DistMatrix<F>& A)¶
Overwrite the uplo triangle of the HPD matrix A with its Cholesky factor.
Note
See HPSDCholesky() for a generalization which also works for semi-definite matrices.
\(LDL^H\) factorization¶
Though the Cholesky factorization is ideal for most HPD matrices, there exist many Hermitian matrices whose eigenvalues are not all positive. The \(LDL^H\) factorization exists as slight relaxation of the Cholesky factorization, i.e., it computes lower-triangular (with unit diagonal) \(L\) and diagonal \(D\) such that \(A = L D L^H\). If \(A\) is found to be numerically singular, then a SingularMatrixException will be thrown.
Warning
The following routines do not pivot, so please use with caution.
- void LDLH(Matrix<F>& A)¶
- void LDLH(DistMatrix<F>& A)¶
Overwrite the strictly lower triangle of \(A\) with the strictly lower portion of \(L\) (\(L\) implicitly has ones on its diagonal) and the diagonal with \(D\).
- void LDLH(Matrix<F>& A, Matrix<F>& d)¶
- void LDLH(DistMatrix<F>& A, DistMatrix<F, MC, STAR>& d)¶
Same as above, but also return the diagonal in the column vector d.
\(LDL^T\) factorization¶
While the \(LDL^H\) factorization targets Hermitian matrices, the \(LDL^T\) factorization targets symmetric matrices. If \(A\) is found to be numerically singular, then a SingularMatrixException will be thrown.
Warning
The following routines do not pivot, so please use with caution.
- void LDLT(Matrix<F>& A)¶
- void LDLT(DistMatrix<F>& A)¶
Overwrite the strictly lower triangle of \(A\) with the strictly lower portion of \(L\) (\(L\) implicitly has ones on its diagonal) and the diagonal with \(D\).
- void LDLT(Matrix<F>& A, Matrix<F>& d)¶
- void LDLT(DistMatrix<F>& A, DistMatrix<F, MC, STAR>& d)¶
Same as above, but also return the diagonal in the vector d.
\(LU\) factorization¶
Given \(A \in \mathbb{F}^{m \times n}\), an LU factorization (without pivoting) finds a unit lower-trapezoidal \(L \in \mathbb{F}^{m \times \mbox{min}(m,n)}\) and upper-trapezoidal \(U \in \mathbb{F}^{\mbox{min}(m,n) \times n}\) such that \(A=LU\). Since \(L\) is required to have its diaganal entries set to one: the upper portion of \(A\) can be overwritten with U, and the strictly lower portion of \(A\) can be overwritten with the strictly lower portion of \(L\). If \(A\) is found to be numerically singular, then a SingularMatrixException will be thrown.
- void LU(Matrix<F>& A)¶
- void LU(DistMatrix<F>& A)¶
Overwrites \(A\) with its LU decomposition.
Since LU factorization without pivoting is known to be unstable for general matrices, it is standard practice to pivot the rows of \(A\) during the factorization (this is called partial pivoting since the columns are not also pivoted). An LU factorization with partial pivoting therefore computes \(P\), \(L\), and \(U\) such that \(PA=LU\), where \(L\) and \(U\) are as described above and \(P\) is a permutation matrix.
- void LU(Matrix<F>& A, Matrix<int>& p)¶
- void LU(DistMatrix<F>& A, DistMatrix<F, VC, STAR>& p)¶
Overwrites the matrix \(A\) with the LU decomposition of \(PA\), where \(P\) is represented by the pivot vector p.
\(LQ\) factorization¶
Given \(A \in \mathbb{F}^{m \times n}\), an LQ factorization typically computes an implicit unitary matrix \(\hat Q \in \mathbb{F}^{n \times n}\) such that \(\hat L \equiv A\hat Q^H\) is lower trapezoidal. One can then form the thin factors \(L \in \mathbb{F}^{m \times \mbox{min}(m,n)}\) and \(Q \in \mathbb{F}^{\mbox{min}(m,n) \times n}\) by setting \(L\) and \(Q\) to first \(\mbox{min}(m,n)\) columns and rows of \(\hat L\) and \(\hat Q\), respectively. Upon completion \(L\) is stored in the lower trapezoid of \(A\) and the Householder reflectors representing \(\hat Q\) are stored within the rows of the strictly upper trapezoid.
- void LQ(Matrix<R>& A)¶
- void LQ(DistMatrix<R>& A)¶
Overwrite the real matrix \(A\) with \(L\) and the Householder reflectors representing \(\hat Q\).
- void LQ(Matrix<Complex<R>>& A, Matrix<Complex<R>>& t)¶
- void LQ(DistMatrix<Complex<R>>& A, DistMatrix<Complex<R>, MD, STAR>& t)¶
Overwrite the complex matrix \(A\) with \(L\) and the Householder reflectors representing \(\hat Q\); unlike the real case, phase information is needed in order to define the (generalized) Householder transformations and is stored in the column vector t.
\(QR\) factorization¶
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.
- void QR(Matrix<R>& A)¶
- void QR(DistMatrix<R>& A)¶
Overwrite the real matrix \(A\) with \(R\) and the Householder reflectors representing \(\hat Q\).
- void QR(Matrix<Complex<R>>& A, Matrix<Complex<R>>& t)¶
- void QR(DistMatrix<Complex<R>>& A, DistMatrix<Complex<R>, MD, STAR>& t)¶
Overwrite the complex matrix \(A\) with \(R\) and the Householder reflectors representing \(\hat Q\); unlike the real case, phase information is needed in order to define the (generalized) Householder transformations and is stored in the column vector t.