LU¶
Driver Routines¶
- template<Number T, DataLocation L, template<typename, typename> typename S1, template<typename, typename> typename S2>
internal::LUSolverDriver<T, L, Matrix<T, L, S2>> lusolve(const Matrix<T, L, S1> &A, const Matrix<T, L, S2> &B)¶Solves the linear system \( \mathbf{A} \; \mathbf{X} = \mathbf{B}\) using LU factorization. ( \( \mathbf{A} = \mathbf{P} \mathbf{L} \mathbf{U}\)). This operation can be done in-place by calling
X = lusolve(A, X)
.
- Parameters:
A – [in] coefficient matrix
B – [in] right-hand side of the system
- Returns:
An abstract object representing the solution of the linear system.
- template<Number T, DataLocation L, template<typename, typename> typename S1, template<typename, typename> typename S2>
internal::LUSolverDriver<T, L, Matrix<T, L, S2>, kInPlace> lusolve(const Matrix<T, L, S1> &A, const Matrix<T, L, S2> &B, InPlace)¶Solves the linear system \( \mathbf{A} \; \mathbf{X} = \mathbf{B}\) using LU factorization ( \( \mathbf{A} = \mathbf{P} \mathbf{L} \mathbf{U}\)). This operation can be done in-place by calling
X = lusolve(A, X)
.Warning
The matrix
A
will be overwritten by the factorL
andU
of the factorization.
- Parameters:
A – [in] coefficient matrix
B – [in] right-hand side of the system
- Returns:
An abstract object representing the solution of the linear system.
- template<Number T, DataLocation L, template<typename, typename> typename S1, template<typename, typename> typename S2>
internal::LUSolverDriver<T, L, Vector<T, L, S2>> lusolve(const Matrix<T, L, S1> &A, const Vector<T, L, S2> &b)¶Solves the linear system \( \mathbf{A} \; \vec{x} = \vec{b}\) using LU factorization ( \( \mathbf{A} = \mathbf{P} \mathbf{L} \mathbf{U}\)). This operation can be done in-place by calling
x = lusolve(A, x)
.
- Parameters:
A – [in] coefficient matrix
b – [in] right-hand side of the system
- Returns:
An abstract object representing the solution of the linear system.
- template<Number T, DataLocation L, template<typename, typename> typename S1, template<typename, typename> typename S2>
internal::LUSolverDriver<T, L, Vector<T, L, S2>> lusolve(const Matrix<T, L, S1> &A, const Vector<T, L, S2> &b, InPlace)¶Solves the linear system \( \mathbf{A} \; \vec{x} = \vec{b}\) using LU factorization ( \( \mathbf{A} = \mathbf{P} \mathbf{L} \mathbf{U}\)). This operation can be done in-place by calling
x = lusolve(A, x)
.Warning
The matrix
A
will be overwritten by the factorL
andU
of the factorization.
- Parameters:
A – [in] coefficient matrix
b – [in] right-hand side of the system
- Returns:
An abstract object representing the solution of the linear system.
- template<Number T, DataLocation L, template<typename, typename> typename S>
LU<T, L> lu(const Matrix<T, L, S> &A)¶Factorize the square matrix
A
as \( \mathbf{A} = \mathbf{PLU}\), whereU
is an upper triangular matrix,L
is a lower triangular matrix with a unit diagonal andP
is a permutation matrix.
- Parameters:
A – [in] input matrix
- template<Number T, DataLocation L, template<typename, typename> typename S>
LU<T, L, kInPlace> lu(const Matrix<T, L, S> &A, InPlace)¶Factorize the square matrix
A
as \( \mathbf{A} = \mathbf{PLU}\), whereU
is an upper triangular matrix,L
is a lower triangular matrix with a unit diagonal andP
is a permutation matrix.Warning
The matrix
A
will be overwritten by the factorL
andU
of the factorization.
- Parameters:
A – [in] input matrix
Factorization¶
-
template<Number T, DataLocation Location, int Flags>
class LU¶ Object representing the LU factorization of a square matrix
A
: \( \mathbf{A} = \mathbf{PLU}\), whereU
is an upper triangular matrix,L
is a lower triangular matrix with a unit diagonal andP
is a permutation matrix.Usage:
Matrix<double> A(n, n); Matrix<double> invA(n, n); Matrix<double> B(n, nrhs); Matrix<double> X(n, nrhs); LU<double> F = lu(A); // Factorize AP = LU X = F.solve(B); // Solve AX = B double cond = F.cond(); // Calculates the reciprocal condition number using L1 norm invA = F.inv(); // Calculate inv(A)
- Template Parameters:
T – Floating-point type
Location – Location of the data.
Flags – Flags for the QR decomposition. Available:
kDefaultFlags
(default flags) andkInPlace
(the factorization is done in-place, overwritingA
).
Public Functions
-
inline LU()¶
Default constructor.
-
template<template<typename, typename> typename S>
inline LU(const Matrix<T, Location, S> &A)¶ Creates a new solver and factorize the matrix
A
.- Parameters:
A – [in] input matrix
-
virtual ~LU() = default¶
Destructor.
-
template<template<typename, typename> typename S>
inline LU &factorize(const Matrix<T, Location, S> &A)¶ Calculates the LU factorization of the matrix
A
.- Parameters:
A – [in] input matrix
-
inline internal::scalar_type_t<T> rcond()¶
- Returns:
The reciprocal of the condition number \( \kappa_1(\mathbf{A}) = \|\mathbf{A}\|_1 \|\mathbf{A}^{-1}\|_1\).
-
template<template<typename, typename> typename S>
inline internal::LUSolve<T, Location, Matrix<T, Location, S>, Flags> solve(Matrix<T, Location, S> &B, char trans = 'N')¶ Solves the linear system \( \text{op}(\mathbf{A}) \; \mathbf{X} = \mathbf{B}\). This operation can be done in-place by calling
X = lu_fact.solve(X)
.- Parameters:
B – [in] right hand side
trans – [in] transpose operation:
N
(no transpose, \( \text{op}(\mathbf{A}) = \mathbf{A} \)),T
(transpose, \( \text{op}(\mathbf{A}) = \mathbf{A}^T \)) orC
(hermitian, complex only, \( \text{op}(\mathbf{A}) = \mathbf{A}^H \))
-
template<template<typename, typename> typename S>
inline internal::LUSolve<T, Location, Vector<T, Location, S>, Flags> solve(Vector<T, Location, S> &b, char trans = 'N')¶ Solves the linear system \( \text{op}(\mathbf{A}) \; \vec{x} = \vec{b}\). This operation can be done in-place by calling
x = lu_fact.solve(x)
.- Parameters:
b – [in] right hand side
trans – [in] transpose operation:
N
(no transpose, \( \text{op}(\mathbf{A}) = \mathbf{A} \)),T
(transpose, \( \text{op}(\mathbf{A}) = \mathbf{A}^T \)) orC
(hermitian, complex only, \( \text{op}(\mathbf{A}) = \mathbf{A}^H \))
-
inline internal::LUInverse<T, Location, Flags> inv()¶
Calculate the inverse of
A
.Warning
To solve the linear system \( \mathbf{A} \; \mathbf{X} = \mathbf{B}\), use call
solve()
instead of explicitly forming the matrix inverse \( \mathbf{A}^{-1}\) and then multiply \( \mathbf{X} = \mathbf{A}^{-1} \mathbf{B}\). This is more efficient and more accurate.
-
inline const mat_type &compressed_mat() const¶
- Returns:
The underlying matrix. The matrix
U
is stored at and above the diagonal, whileL
is stored below the diagonal. The unit diagonal ofL
is not stored.
-
inline TriangularMap<T, kLowerUnit, Location> L() const¶
- Returns:
The lower triangular matrix
L
from the LU factorization. The unit diagonal ofL
is not stored.
-
inline TriangularMap<T, kUpper, Location> U() const¶
- Returns:
The upper triangular matrix
U
from the LU factorization.