Updateable Cholesky decomposition: A=LDL'. More...
#include <vnl_ldl_cholesky.h>
Public Types | |
| enum | Operation { quiet, verbose, estimate_condition } |
| Modes of computation. See constructor for details. More... | |
Public Member Functions | |
| vnl_ldl_cholesky (vnl_matrix< double > const &M, Operation mode=verbose) | |
| Make cholesky decomposition of M optionally computing the reciprocal condition number. | |
| ~vnl_ldl_cholesky () | |
| vnl_vector< double > | solve (vnl_vector< double > const &b) const |
| Solve LS problem M x = b. | |
| void | solve (vnl_vector< double > const &b, vnl_vector< double > *x) const |
| Solve LS problem M x = b. | |
| void | solve_lx (vnl_vector< double > &y) |
| Solve equation of form Lx=y (in-place). | |
| double | determinant () const |
| Compute determinant. | |
| void | rank1_update (const vnl_vector< double > &v) |
| Compute rank-1 update, ie the decomposition of (M+v.v'). | |
| void | update (const vnl_matrix< double > &W) |
| Multi-rank update, ie the decomposition of (M+W.W'). | |
| vnl_matrix< double > | inverse () const |
| Compute inverse. Not efficient. | |
| const vnl_matrix< double > & | lower_triangle () const |
| Return lower-triangular factor. | |
| vnl_matrix< double > | upper_triangle () const |
| Return upper-triangular factor. | |
| const vnl_vector< double > & | diagonal () const |
| Return elements of diagonal matrix D in LDL'. | |
| double | xt_m_inv_x (const vnl_vector< double > &x) const |
| Efficient computation of x' * inv(M) * x. | |
| double | xt_m_x (const vnl_vector< double > &x) const |
| Efficient computation of x' * M * x. | |
| int | rank_deficiency () const |
| A Success/failure flag. | |
| double | rcond () const |
| Return reciprocal condition number (smallest/largest singular values). | |
| vnl_vector< double > & | nullvector () |
| Return computed nullvector. | |
| vnl_vector< double > const & | nullvector () const |
Protected Attributes | |
| vnl_matrix< double > | L_ |
| Lower triangular matrix. | |
| vnl_vector< double > | d_ |
| Elements of diagonal matrix. | |
| double | rcond_ |
| 1/(condition number). | |
| long | num_dims_rank_def_ |
| vnl_vector< double > | nullvector_ |
Private Member Functions | |
| vnl_ldl_cholesky (vnl_ldl_cholesky const &that) | |
| Copy constructor - privatised to avoid it being used. | |
| vnl_ldl_cholesky & | operator= (vnl_ldl_cholesky const &that) |
| Assignment operator - privatised to avoid it being used. | |
| void | inplace_solve (double *x) const |
| Solve Mx=b, overwriting input vector with the solution. | |
Updateable Cholesky decomposition: A=LDL'.
A class to hold the Cholesky decomposition of a positive definite symmetric matrix, using the form A=LDL', where L is lower triangular with ones on the leading diagonal, and D is a diagonal matrix. This differs from vnl_cholesky, which decomposes as A=LL', without a constraint on the diagonal. The advantage of the LDL' form is that it can be efficiently updated. Thus given the decomposition of A, we can compute the decomposition of (A+v.v') in O(n^2) time.
This can be used to solve linear systems, compute determinants and inverses.
To check that the decomposition can be used safely for solving a linear equation it is wise to construct with mode==estimate_condition and check that rcond()>sqrt(machine precision). If this is not the case it might be a good idea to use vnl_svd instead.
Definition at line 31 of file vnl_ldl_cholesky.h.
Modes of computation. See constructor for details.
Definition at line 35 of file vnl_ldl_cholesky.h.
| vnl_ldl_cholesky::vnl_ldl_cholesky | ( | vnl_matrix< double > const & | M, |
| Operation | mode = verbose |
||
| ) |
Make cholesky decomposition of M optionally computing the reciprocal condition number.
Cholesky decomposition.
Make cholesky decomposition of M optionally computing the reciprocal condition number. If mode is estimate_condition, the condition number and an approximate nullspace are estimated, at a cost of a factor of (1 + 18/n). Here's a table of 1 + 18/n:
n: 3 5 10 50 100 500 1000 slowdown: 7.0 4.6 2.8 1.4 1.18 1.04 1.02
Sqrt of elements of diagonal matrix.
Definition at line 28 of file vnl_ldl_cholesky.cxx.
| vnl_ldl_cholesky::~vnl_ldl_cholesky | ( | ) | [inline] |
Definition at line 43 of file vnl_ldl_cholesky.h.
| vnl_ldl_cholesky::vnl_ldl_cholesky | ( | vnl_ldl_cholesky const & | that | ) | [private] |
Copy constructor - privatised to avoid it being used.
| double vnl_ldl_cholesky::determinant | ( | ) | const |
Compute determinant.
Definition at line 180 of file vnl_ldl_cholesky.cxx.
| const vnl_vector<double>& vnl_ldl_cholesky::diagonal | ( | ) | const [inline] |
Return elements of diagonal matrix D in LDL'.
Definition at line 80 of file vnl_ldl_cholesky.h.
| void vnl_ldl_cholesky::inplace_solve | ( | double * | x | ) | const [private] |
Solve Mx=b, overwriting input vector with the solution.
x points to beginning of an n-element vector containing b On exit, x[i] filled with solution vector.
Definition at line 101 of file vnl_ldl_cholesky.cxx.
| vnl_matrix< double > vnl_ldl_cholesky::inverse | ( | ) | const |
Compute inverse. Not efficient.
Note that you rarely need the inverse - backsubstitution is faster and less prone to rounding errors.
Definition at line 252 of file vnl_ldl_cholesky.cxx.
| const vnl_matrix<double>& vnl_ldl_cholesky::lower_triangle | ( | ) | const [inline] |
Return lower-triangular factor.
Definition at line 74 of file vnl_ldl_cholesky.h.
| vnl_vector<double>& vnl_ldl_cholesky::nullvector | ( | ) | [inline] |
Return computed nullvector.
Not calculated unless Operation mode at construction was estimate_condition.
Definition at line 101 of file vnl_ldl_cholesky.h.
| vnl_vector<double> const& vnl_ldl_cholesky::nullvector | ( | ) | const [inline] |
Definition at line 102 of file vnl_ldl_cholesky.h.
| vnl_ldl_cholesky& vnl_ldl_cholesky::operator= | ( | vnl_ldl_cholesky const & | that | ) | [private] |
Assignment operator - privatised to avoid it being used.
| void vnl_ldl_cholesky::rank1_update | ( | const vnl_vector< double > & | v | ) |
Compute rank-1 update, ie the decomposition of (M+v.v').
If the initial state is the decomposition of M, then L and D are updated so that on exit LDL'=M+v.v'
If the initial state is the decomposition of M, then L and D are updated so that on exit LDL'=M+v.v'
Uses the algorithm given by Davis and Hager in "Multiple-Rank Modifications of a Sparse Cholesky Factorization",2001.
Definition at line 194 of file vnl_ldl_cholesky.cxx.
| int vnl_ldl_cholesky::rank_deficiency | ( | ) | const [inline] |
A Success/failure flag.
Definition at line 91 of file vnl_ldl_cholesky.h.
| double vnl_ldl_cholesky::rcond | ( | ) | const [inline] |
Return reciprocal condition number (smallest/largest singular values).
As long as rcond()>sqrt(precision) the decomposition can be used for solving equations safely. Not calculated unless Operation mode at construction was estimate_condition.
Definition at line 97 of file vnl_ldl_cholesky.h.
| vnl_vector< double > vnl_ldl_cholesky::solve | ( | vnl_vector< double > const & | b | ) | const |
Solve LS problem M x = b.
Solve least squares problem M x = b.
Definition at line 170 of file vnl_ldl_cholesky.cxx.
| void vnl_ldl_cholesky::solve | ( | vnl_vector< double > const & | b, |
| vnl_vector< double > * | xp | ||
| ) | const |
Solve LS problem M x = b.
Solve least squares problem M x = b.
The right-hand-side vcl_vector x may be b, which will give a fractional increase in speed.
Definition at line 161 of file vnl_ldl_cholesky.cxx.
| void vnl_ldl_cholesky::solve_lx | ( | vnl_vector< double > & | x | ) |
Solve equation of form Lx=y (in-place).
Solve Lx=y (in-place).
x is overwritten with solution
Definition at line 91 of file vnl_ldl_cholesky.cxx.
| void vnl_ldl_cholesky::update | ( | const vnl_matrix< double > & | W0 | ) |
Multi-rank update, ie the decomposition of (M+W.W').
If the initial state is the decomposition of M, then L and D are updated so that on exit LDL'=M+W.W'
Definition at line 219 of file vnl_ldl_cholesky.cxx.
| vnl_matrix<double> vnl_ldl_cholesky::upper_triangle | ( | ) | const [inline] |
Return upper-triangular factor.
Definition at line 77 of file vnl_ldl_cholesky.h.
| double vnl_ldl_cholesky::xt_m_inv_x | ( | const vnl_vector< double > & | x | ) | const |
Efficient computation of x' * inv(M) * x.
Useful when M is a covariance matrix!
Useful when M is a covariance matrix! Solves Ly=x for y, then returns sum y[i]*y[i]/d[i]
Definition at line 124 of file vnl_ldl_cholesky.cxx.
| double vnl_ldl_cholesky::xt_m_x | ( | const vnl_vector< double > & | x | ) | const |
Efficient computation of x' * M * x.
Twice as fast as explicitly computing x' * M * x
Definition at line 141 of file vnl_ldl_cholesky.cxx.
vnl_vector<double> vnl_ldl_cholesky::d_ [protected] |
Elements of diagonal matrix.
Definition at line 111 of file vnl_ldl_cholesky.h.
vnl_matrix<double> vnl_ldl_cholesky::L_ [protected] |
Lower triangular matrix.
Definition at line 108 of file vnl_ldl_cholesky.h.
vnl_vector<double> vnl_ldl_cholesky::nullvector_ [protected] |
Definition at line 116 of file vnl_ldl_cholesky.h.
long vnl_ldl_cholesky::num_dims_rank_def_ [protected] |
Definition at line 115 of file vnl_ldl_cholesky.h.
double vnl_ldl_cholesky::rcond_ [protected] |
1/(condition number).
Definition at line 114 of file vnl_ldl_cholesky.h.
1.7.5.1