Interval LU decomposition
Main author: Damien Massé
This page relates to the LU decomposition of an interval matrix with complete pivoting, and related features.
-
class IntvFullPivLU
Full pivot LU decomposition for a matrix of intervals, based on Eigen decomposition. The decomposition is of the form \(\mathbf{M} = \mathbf{P}^{-1} [\mathbf{L}][\mathbf{U}] \mathbf{Q}^{-1}\) where \(\mathbf{P}\) and \(\mathbf{Q}\) are permutation matrices, and \([\mathbf{L}]\) and \([\mathbf{U}]\) are lower and upper interval matrices (i.e. \([\mathbf{L}]\)’s diagonal is 1).
Public Functions
-
explicit IntvFullPivLU(const Matrix &M)
Constructor from
Matrixof double.- Parameters:
M – the matrix for which the decomposition is computed
-
explicit IntvFullPivLU(const IntervalMatrix &M)
Constructor from
Matrixof intervals. Eigen decomposition is done onM.mid().- Parameters:
M – the matrix of intervals.
-
BoolInterval is_injective() const
Check if the matrix is injective, i.e. its rank is equal to its number of rows.
- Returns:
-
BoolInterval is_invertible() const
Check if the initial matrix is invertible i.e. it is square and full rank.
- Returns:
-
BoolInterval is_surjective() const
Check if the matrix is surjective i.e. its rank is equal to its number of cols.
- Returns:
-
Interval determinant() const
Return an interval enclosing the determinant.
- Pre:
The matrix must be square.
- Returns:
the product of the diagonal elements of \([\mathbf{U}]\)
-
Interval rank() const
Return an interval enclosing the rank. Quite precise for square matrix (number of diagonal elements of \([\mathbf{U}]\) not containing 0). Less precise for non-square matrices, where each row/column outside the top-left part of \([\mathbf{U}]\) can change the rank. However, if no diagonal element contains 0, the return is unambiguous.
- Returns:
an interval enclosing the possible ranks.
-
Interval dimension_of_kernel() const
Approximation of the size of the kernel space, based on the result of
rank()(number of cols-rank()). As such, this is not the exact size of the kernel space as built bykernel().- Returns:
an interval enclosing the possible dimensions.
-
IntervalMatrix kernel() const
Overapproximation of the kernel space as a matrix of column vectors. Any vector \(\mathbf{V}\) which is not a linear combination of the column vectors is guaranteed to be outside the kernel (i.e., \(\mathbf{M}\mathbf{V} \neq 0\)).
- Returns:
a matrix of column vectors, may be empty
-
IntervalMatrix cokernel() const
Overapproximation of the left-null (“cokernel”) space as a matrix of row vectors. Any vector \(\mathbf{V}\) which is not a linear combination of the row vectors is guaranteed to be outside the kernel (i.e., \(\mathbf{V}\mathbf{M} \neq 0\)).
- Returns:
a matrix of row vectors, may be empty
-
template<typename Derived>
inline Derived image(const Eigen::MatrixBase<Derived> &M) const “Underapproximation” of the column space of the matrix, i.e. return a set of independant columns of the original matrix which is possibly maximal. As for Eigen, you must provide the original matrix used for the decomposition.
- Parameters:
M – the matrix used to build the decomposition.
- Pre:
\(\mathbf{M}\) is the matrix used to build the decomposition
- Returns:
a matrix of columns of \(\mathbf{M}\), or one vector of 0 if the rank of \(\mathbf{M}\) may be 0.
-
template<typename Derived>
inline Derived coimage(const Eigen::MatrixBase<Derived> &M) const “Underapproximation” of the row space of the matrix, i.e. return a set of independant rows of the original matrix which is possibly maximal. As for Eigen, you must provide the original matrix used for the decomposition.
- Parameters:
M – the matrix used to build the decomposition.
- Pre:
M is the matrix used to build the decomposition
- Returns:
a matrix of rows of \(\mathbf{M}\), or one row of 0 if the rank of \(\mathbf{M}\) may be 0.
-
IntervalMatrix solve(const IntervalMatrix &rhs) const
Equation solving \(\mathbf{M}\mathbf{X}=\mathbf{rhs}\).
Precisely look for solutions where the only non-zero values are those on non-zero pivots. If the matrix is full-rank and surjective (cols >= rows), it gives an overapproximation of the solutions (for each possible values of rhs).
If the matrix is full-rank and injective (rows >= cols), it returns empty if no solution is possible, and a possible overapproximation of the solutions otherwise (but there may still be no solution).
If the matrix is not full-rank,
empty means that no solution is possible with the initial precondition (non-zero values for non-zero pivots) if the goal is to prove the absence of solution, see
solvewith a bounding box.non empty presents the possible solutions found.
- Parameters:
rhs – right-hand side of the equation
- Returns:
a potential solution of the equation \(\mathbf{M}\mathbf{X}=\mathbf{rhs}\)
-
void solve(const IntervalMatrix &rhs, IntervalMatrix &B) const
Equation solving \(\mathbf{M}\mathbf{X}=\mathbf{rhs}\) with bounding matrix \([\mathbf{B}]\) for the solution, i.e. contraction of the matrix on the solutions on \(\mathbf{M}\mathbf{X}\).
Especially useful where the matrix is not full-rank, where solve without bounding box may return empty even if solutions exist.
If the matrix is full-rank and surjective (cols >= rows), it gives an overapproximation of the solutions (for each possible values of rhs).
If the matrix is full-rank and injective (rows >= cols), it returns empty if no solution is possible, and a possible overapproximation of the solutions otherwise (but there may still be no solution).
If the matrix is not full-rank,
empty means that no solution is possible inside the bounding matrix.
non empty presents the possible solutions found.
- Parameters:
rhs – right-hand side of the equation
B – bounding-box of the left hand-side, contracted so that \(\mathbf{M}\mathbf{B}=\mathbf{rhs}\)
-
IntervalMatrix reconstructed_matrix() const
Rebuilding of the matrix, i.e. compute \(\mathbf{P}^{-1}[\mathbf{L}][\mathbf{U}]\mathbf{Q}^{-1}\).
Can be used to evaluate the precision of the decomposition.
- Returns:
the reconstructed matrix
-
double max_pivot() const
Maximum magnitude of the diagonal elements of \([\mathbf{U}]\).
- Returns:
the maximum
-
inline const Eigen::FullPivLU<Matrix>::PermutationPType &permutation_P() const
The permutation \(\mathbf{P}\) in the decomposition \(\mathbf{P}^{-1}\mathbf{L}\mathbf{U}\mathbf{Q}^{-1}\).
- Returns:
the permutation \(\mathbf{P}\), as defined by Eigen
-
inline const Eigen::FullPivLU<Matrix>::PermutationQType &permutation_Q() const
The permutation \(\mathbf{Q}\) in the decomposition \(\mathbf{P}^{-1}\mathbf{L}\mathbf{U}\mathbf{Q}^{-1}\).
- Returns:
the permutation \(\mathbf{Q}\), as defined by Eigen
-
inline const Eigen::FullPivLU<Matrix> &eigen_LU() const
The Eigen decomposition of
M.mid()- Returns:
the Eigen decomposition
-
inline const Row &transformation() const
Return the column-wise transformation done on
M.mid()before the Eigen LU decomposition, as a Row.- Returns:
the transformation
-
inline const IntervalMatrix &matrix_LU() const
Returns the matrix storing \([\mathbf{L}]\) and \([\mathbf{U}]\) (i.e. \([\mathbf{L}]\) for strictly lower part, \([\mathbf{U}]\) for upper part).
- Returns:
the interval matrix
-
explicit IntvFullPivLU(const Matrix &M)