Matrix Decomposition Methods
The following methods are available for matrix decompositions and solution of simultaneous linear equations.
1. T lud (valarray<size_t>& ri);
2. void lubksb (const valarray<size_t>& ri, const valarray<T>& v, valarray<T>& s) const;
3. void lumpove (const matrix<T>& ludm, const valarray<size_t>& ri, const valarray<T>& v, valarray<T>& s) const;
4. void svd (matrix<T>& z, valarray<T>& w);
5. void svbksb (const matrix<T>& z, const valarray<T>& w, const valarray<T>& v, valarray<T>& s) const;
6. void qrd (matrix<T>& r);
7. void qrbksb (const matrix<T>& r, const valarray<T>& v, valarray<T>& s) const;
8. void chold ();
9. void cholbksb (const valarray<T>& v, valarray<T>& s) const;
10. void inv ();
11. void inv2 ();
12. void solve (const valarray<T>& v, valarray<T>& s) const;
13. void solve (const matrix<T>& v, matrix<T>& s) const;
- Performs LU decomposition on a square matrix, i.e., the matrix A is decomposed as,
A = L . U
Where L is a lower triangular matrix and U is an upper triangular matrix. On decomposition, the matrix A is replaced by both L and U. The new changed row index is returned in the valarray ri for performing the back substitution later. The method throws invalid_argument exception for non-square matrix if range checking is enabled. It can also throw runtime_error exception if the matrix is found singular during decomposition. Returns the determinant of the matrix.
- Performs forward and back substitution on a LU decomposed matrix to find out the solution vector for a linear set of equations. Here, v is the right hand side input vector and ri is the row index as returned by lud method. The solution vector will be returned in the vector s.
- The lumpove method can be used to improve an imperfect solutions found by LU decomposition method earlier. This method must be applied to the original matrix. Here, ludm is the LU decomposition of the original matrix, v is the right hand side vector and ri is the row index. The improved solution will be returned in the vector s.
- Performs Singular Value decomposition on a matrix, i.e., any M x N matrix A, where M >= N, is decomposed as,
A = U . W . ZT
Where U is an M x N column-orthogonal matrix, W is a diagonal matrix, and Z is an N x N orthogonal matrix. After successful decomposition, the original matrix A is replaced by U, the diagonal matrix is returned in the vector w, and the orthogonal matrix in z. This method may throw exception runtime_error if there is no convergence after 30 iterations.
- The svbksb method performs back substitution on a SV decomposed matrix to find out the solution vector for a linear set of equations. Here the matrix z and vector w are as returned by the svd method, v is the right hand side vector, and s is the solution vector to be returned by this method.
- The QR decomposition method decomposes an M-by-N matrix A, where M >= N, into an M-by-N orthogonal matrix Q and an N-by-N upper triangular matrix R, i.e.,
A = Q . R
After successful decomposition, the original matrix A is replaced by the orthogonal matrix Q, and the upper triangular matrix R is returned in r.
- Finds the solution vector for a QR decomposed matrix. Here r is the upper triangular matrix returned by the qrd method, v is the right hand side vector, and s is the reference to a valarray class object to receive the solution vector.
- The Cholesky decomposition method chold decomposes a symmetric, positive definite matrix A, into a lower triangular matrix L so that
A = L . LT
It will throw runtime_error exception if the matrix is not symmetric and positive definite.
- The cholbksb method performs back substitution on a Cholesky decomposed matrix to find out the solution vector for a set of linear equations. Here v is the right hand side vector, and s is the reference to valarray object to receive the solution vector.
- Inverts the matrix. The matrix is inverted in-place using minimum memory, i.e., the original matrix is replaced by its inverse.
- Inverts the matrix. The matrix is inverted using LU decomposition method. The original matrix is replaced by its inverse. This method uses more memory than inv and is also slightly slower.
- Solves a set of linear equations using LU decomposition method. Here v is the right hand side vector, and s is the reference to valarray object to receive the solution vector.
- Solves a set of linear equations using LU decomposition method. The matrix v can contain more than one the right hand side vector as column. For each right hand side vector, the solution vector is returned in the corresponding column of the matrix s.
Examples
typedef matrix<double> Matrix;
typedef valarray<double> Vector;
Matrix A(5,5), Z(5,5);
Vector V(5), S(5), W(5);
valarray<size_t> ri;
cin >> A >> V; // Reads the matrix and vector from stdin
Matrix B = A;
A.lud( ri);
A.lubksb( ri, V, S);
A = B;
A.svd( Z, W);
A.svbksb( Z, W, V, S);
A = B;
A.qrd( Z);
A.qrbksb( Z, V, S);
A = B[gmslice(LTRIANG)];
A[gmslice(LTRIANG)] = ~A; // A is now a symmetric matrix
A.chold();
A.cholbksb( V, S);
A = B;
A.inv(); // Inverts the matrix
A = B;
A.solve( V, S);