QRMatrix (i.e. QR decomposition, QR factorisation or orthogonal-triangular decomposition) decomposes a scalar/complex matrix A
into the following matrix product:
More...
Public Types | |
enum | outputTypes : uint8_t { FULL_R = 1, FULL_QR = 2, REDUCED_R = 3 } |
Options for the output matrix forms of QRMatrix. More... | |
enum | storeMethods : uint8_t { IN_PLACE = 1, OUT_OF_PLACE = 2 } |
Options where to store R. More... | |
enum | colPivoting : bool { FALSE = false, TRUE = true } |
Options for the computation of column pivoting. More... | |
typedef MatrixType::cmptType | cmptType |
typedef SquareMatrix< cmptType > | SMatrix |
typedef RectangularMatrix< cmptType > | RMatrix |
Public Member Functions | |
QRMatrix () | |
Construct null. More... | |
QRMatrix (const outputTypes outputType, const storeMethods storeMethod=storeMethods::OUT_OF_PLACE, const colPivoting colPivot=colPivoting::FALSE) | |
Construct QRMatrix without performing the decomposition. More... | |
QRMatrix (MatrixType &A, const outputTypes outputType, const storeMethods storeMethod=storeMethods::OUT_OF_PLACE, const colPivoting colPivot=colPivoting::FALSE) | |
Construct QRMatrix and perform the QR decomposition. More... | |
QRMatrix (const MatrixType &A, const outputTypes outputType, const colPivoting colPivot=colPivoting::FALSE) | |
Construct QRMatrix and perform the QR decomposition. More... | |
const MatrixType & | Q () const |
Return the unitary similarity matrix. More... | |
const MatrixType & | R () const |
Return the upper triangular matrix. More... | |
const labelList & | orderP () const |
Return the permutation order (P) as a list. More... | |
SMatrix | P () const |
Create and return the permutation matrix. More... | |
void | decompose (MatrixType &A) |
Compute QR decomposition according to constructor settings. More... | |
void | decompose (const MatrixType &A) |
Compute QR decomposition according to constructor settings. More... | |
void | solve (List< cmptType > &x, const UList< cmptType > &source) const |
template<class Addr > | |
void | solve (List< cmptType > &x, const IndirectListBase< cmptType, Addr > &source) const |
tmp< Field< cmptType > > | solve (const UList< cmptType > &source) const |
template<class Addr > | |
tmp< Field< cmptType > > | solve (const IndirectListBase< cmptType, Addr > &source) const |
SMatrix | inv () const |
Return the inverse of (Q*R), so that solving x = (Q*R).inv()*source. More... | |
RMatrix | backSubstitution (const SMatrix &A, const RMatrix &b) |
Solve a row-echelon-form linear system starting from the bottom. More... | |
template<class Addr > | |
Foam::tmp< Foam::Field< typename MatrixType::cmptType > > | solve (const IndirectListBase< cmptType, Addr > &source) const |
Protected Member Functions | |
RMatrix | householderReflector (RMatrix u) |
Compute Householder reflector on a given matrix column, u. More... | |
void | applyLeftReflector (MatrixType &A, const RMatrix &reflector, const label k=0, const label k1=0) |
Apply (in-place) Householder reflectors from the left side: u*A. More... | |
void | applyRightReflector (MatrixType &A, const RMatrix &reflector, const label k=0) |
Apply (in-place) Householder reflectors from the right side: (u*A)*u. More... | |
QRMatrix (i.e. QR decomposition, QR factorisation or orthogonal-triangular decomposition) decomposes a scalar/complex matrix A
into the following matrix product:
A = Q*R,
where Q
is a unitary similarity matrix, R
is an upper triangular matrix.
Reference:
QR decomposition: mathworld.wolfram.com/QRDecomposition.html (Retrieved:17-06-19) w.wiki/52X (Retrieved: 17-06-19) QR decomposition with Householder reflector (tag:M): Monahan, J. F. (2011). Numerical methods of statistics. Cambridge: Cambridge University Press. DOI:10.1017/CBO9780511977176 QR decomposition with column pivoting (tag:QSB): Quintana-Ortí, G., Sun, X., & Bischof, C. H. (1998). A BLAS-3 version of the QR factorization with column pivoting. SIAM Journal on Scientific Computing, 19(5), 1486-1494. DOI:10.1137/S1064827595296732 Moore-Penrose inverse algorithm (tags:KP; KPP): Katsikis, V. N., & Pappas, D. (2008). Fast computing of the Moore-Penrose inverse matrix. Electronic Journal of Linear Algebra, 17(1), 637-650. DOI:10.13001/1081-3810.1287 Katsikis, V. N., Pappas, D., & Petralias, A. (2011). An improved method for the computation of the Moore–Penrose inverse matrix. Applied Mathematics and Computation, 217(23), 9828-9834. DOI:10.1016/j.amc.2011.04.080 Tolerance for the Moore-Penrose inverse algorithm (tag:TA): Toutounian, F., & Ataei, A. (2009). A new method for computing Moore–Penrose inverse matrices. Journal of Computational and applied Mathematics, 228(1), 412-417. DOI:10.1016/j.cam.2008.10.008 Operands of QR decomposition: mathworld.wolfram.com/UnitaryMatrix.html (Retrieved:13-06-19) mathworld.wolfram.com/UpperTriangularMatrix.html (Retrieved:16-06-19) mathworld.wolfram.com/PermutationMatrix.html (Retrieved:20-06-19) Back substitution: mathworld.wolfram.com/GaussianElimination.html (Retrieved:15-06-19)
A
can be a SquareMatrix<Type>
or RectangularMatrix<Type>
Output types:
Q
is always of the type of the matrix A
R
is always of the type of the matrix A
Options for the output forms of QRMatrix
(for an (m-by-n) input matrix A
with k = min(m, n)):
R
(m-by-n)R
and Q
(m-by-m)R
(k-by-n)Options where to store R:
R
R
Options for the computation of column pivoting:
Direct solution of linear systems A x = b is possible by solve() alongside the following limitations:
A
= a scalar square matrixNotes
R
is not positive diagonal R
.A
after its nth row.R
if m > n for an input m-by-n matrix A
.Definition at line 152 of file QRMatrix.H.
typedef MatrixType::cmptType cmptType |
Definition at line 157 of file QRMatrix.H.
typedef SquareMatrix<cmptType> SMatrix |
Definition at line 158 of file QRMatrix.H.
typedef RectangularMatrix<cmptType> RMatrix |
Definition at line 159 of file QRMatrix.H.
enum outputTypes : uint8_t |
Options for the output matrix forms of QRMatrix.
Enumerator | |
---|---|
FULL_R | computes only |
FULL_QR | computes both |
REDUCED_R | computes only reduced |
Definition at line 162 of file QRMatrix.H.
enum storeMethods : uint8_t |
Options where to store R.
Enumerator | |
---|---|
IN_PLACE | replaces input matrix content with |
OUT_OF_PLACE | creates new object of |
Definition at line 170 of file QRMatrix.H.
enum colPivoting : bool |
Options for the computation of column pivoting.
Enumerator | |
---|---|
FALSE | switches off column pivoting |
TRUE | switches on column pivoting |
Definition at line 177 of file QRMatrix.H.
QRMatrix | ( | ) |
Construct null.
Definition at line 186 of file QRMatrix.C.
QRMatrix | ( | const outputTypes | outputType, |
const storeMethods | storeMethod = storeMethods::OUT_OF_PLACE , |
||
const colPivoting | colPivot = colPivoting::FALSE |
||
) |
Construct QRMatrix without performing the decomposition.
Definition at line 196 of file QRMatrix.C.
QRMatrix | ( | MatrixType & | A, |
const outputTypes | outputType, | ||
const storeMethods | storeMethod = storeMethods::OUT_OF_PLACE , |
||
const colPivoting | colPivot = colPivoting::FALSE |
||
) |
Construct QRMatrix and perform the QR decomposition.
Definition at line 213 of file QRMatrix.C.
References A.
QRMatrix | ( | const MatrixType & | A, |
const outputTypes | outputType, | ||
const colPivoting | colPivot = colPivoting::FALSE |
||
) |
Construct QRMatrix and perform the QR decomposition.
Definition at line 233 of file QRMatrix.C.
References A.
|
inlineprotected |
Compute Householder reflector on a given matrix column, u.
(M:Eq.6.814)
Definition at line 52 of file QRMatrixI.H.
References Foam::abort(), Matrix< Form, Type >::columnNorm(), Foam::exit(), Foam::FatalError, FatalErrorInFunction, Foam::mag(), and Matrix< Form, Type >::n().
|
inlineprotected |
Apply (in-place) Householder reflectors from the left side: u*A.
Definition at line 97 of file QRMatrixI.H.
References A, Foam::Detail::conj(), k, Matrix< Form, Type >::m(), Foam::sum(), and Foam::Zero.
|
inlineprotected |
Apply (in-place) Householder reflectors from the right side: (u*A)*u.
Definition at line 126 of file QRMatrixI.H.
References A, Foam::Detail::conj(), k, Matrix< Form, Type >::m(), Foam::sum(), and Foam::Zero.
|
inline |
Return the unitary similarity matrix.
Includes implicit round-to-zero as mutable operation
Definition at line 155 of file QRMatrixI.H.
Referenced by Foam::pinv().
|
inline |
Return the upper triangular matrix.
Definition at line 164 of file QRMatrixI.H.
Referenced by Foam::pinv().
|
inline |
Return the permutation order (P) as a list.
Definition at line 171 of file QRMatrixI.H.
|
inline |
Create and return the permutation matrix.
Definition at line 179 of file QRMatrixI.H.
References forAll, and Foam::Zero.
Referenced by Foam::pinv().
void decompose | ( | MatrixType & | A | ) |
Compute QR decomposition according to constructor settings.
Definition at line 254 of file QRMatrix.C.
References A, Foam::I, Foam::nl, and WarningInFunction.
void decompose | ( | const MatrixType & | A | ) |
Compute QR decomposition according to constructor settings.
Definition at line 313 of file QRMatrix.C.
References A, Foam::I, Foam::nl, and WarningInFunction.
Solve the linear system with the given source and return the solution in the argument x
Definition at line 358 of file QRMatrix.C.
References x.
void solve | ( | List< cmptType > & | x, |
const IndirectListBase< cmptType, Addr > & | source | ||
) | const |
Solve the linear system with the given source and return the solution in the argument x
Definition at line 370 of file QRMatrix.C.
References x.
Foam::tmp< Foam::Field< typename MatrixType::cmptType > > solve | ( | const UList< cmptType > & | source | ) | const |
Solve the linear system with the given source and return the solution
Definition at line 382 of file QRMatrix.C.
tmp<Field<cmptType> > solve | ( | const IndirectListBase< cmptType, Addr > & | source | ) | const |
Solve the linear system with the given source and return the solution
Foam::QRMatrix< MatrixType >::SMatrix inv | ( | ) | const |
Return the inverse of (Q*R), so that solving x = (Q*R).inv()*source.
Definition at line 412 of file QRMatrix.C.
References Foam::inv(), Matrix< Form, Type >::m(), and x.
Foam::RectangularMatrix< typename MatrixType::cmptType > backSubstitution | ( | const SMatrix & | A, |
const RMatrix & | b | ||
) |
Solve a row-echelon-form linear system starting from the bottom.
Back substitution: Ax = b where: A = Non-singular upper-triangular square matrix (m-by-m) x = Solution (m-by-p) b = Source (m-by-p)
Definition at line 436 of file QRMatrix.C.
References A, Foam::abort(), Foam::constant::atomic::alpha, Foam::diag(), Foam::FatalError, FatalErrorInFunction, k, Matrix< Form, Type >::m(), Foam::mag(), Matrix< Form, Type >::n(), Foam::nl, Foam::tab, WarningInFunction, and Foam::Zero.
Referenced by Foam::pinv().
Foam::tmp<Foam::Field<typename MatrixType::cmptType> > solve | ( | const IndirectListBase< cmptType, Addr > & | source | ) | const |
Definition at line 398 of file QRMatrix.C.