template<typename _MatrixType, int _UpLo>
class Eigen::LDLT< _MatrixType, _UpLo >
Robust Cholesky decomposition of a matrix with pivoting.
- Parameters
-
MatrixType | the type of the matrix of which to compute the LDL^T Cholesky decomposition |
Perform a robust Cholesky decomposition of a positive semidefinite or negative semidefinite matrix
such that
, where P is a permutation matrix, L is lower triangular with a unit diagonal and D is a diagonal matrix.
The decomposition uses pivoting to ensure stability, so that L will have zeros in the bottom right rank(A) - n submatrix. Avoiding the square root on D also stabilizes the computation.
Remember that Cholesky decompositions are not rank-revealing. Also, do not use a Cholesky decomposition to determine whether a system of equations has a solution.
- See Also
- MatrixBase::ldlt(), class LLT
const internal::solve_retval<LDLT, Rhs> solve |
( |
const MatrixBase< Rhs > & |
b | ) |
const |
|
inline |
- Returns
- a solution x of
using the current decomposition of A.
This function also supports in-place solves using the syntax x = decompositionObject.solve(x)
.
This method just tries to find as good a solution as possible. If you want to check whether a solution exists or if it is accurate, just call this function to get a result and then compute the error of this result, or use MatrixBase::isApprox() directly, for instance like this:
bool a_solution_exists = (A*result).isApprox(b, precision);
This method avoids dividing by zero, so that the non-existence of a solution doesn't by itself mean that you'll get inf
or nan
values.
More precisely, this method solves
using the decomposition
by solving the systems
,
,
,
and
in succession. If the matrix
is singular, then
will also be singular (all the other matrices are invertible). In that case, the least-square solution of
is computed. This does not mean that this function computes the least-square solution of
is
is singular.
- See Also
- MatrixBase::ldlt()