#include <rw/flsch.h> // FloatLeastSqCh, float Cholesky method #include <rw/dlsch.h> // DoubleLeastSqCh, double Cholesky method #include <rw/clsch.h> // DComplexLeastSqCh, complex Cholesky
// method #include <rw/flsqr.h> // FloatLeastSqQR, float QR method #include <rw/dlsqr.h> // DoubleLeastSqQR, double QR method #include <rw/clsqr.h> // DComplexLeastSqQR, complex QR method #include <rw/flssv.h> // FloatLeastSqSV, float SV method #include <rw/dlssv.h> // DoubleLeastSqSV, double SV method #include <rw/clssv.h> // DComplexLeastSqSV, complex SV method DoubleLeastSqSV SVD(A); // A is an RWGenMat<double> RWMathVec<double> x = SVD.solve(b); // b is an RWMathVec<double>
A linear system of equations has the form Ax=b, where A contains more rows than columns and generally has no exact solution; it is overdetermined. The best that we can do is to find a solution vector x such that the norm of the residual vector, Ax-b, is made as small as possible. The vector x is then a solution in the least squares sense, and this sort of problem is called a least squares problem.
As shown in the Synopsis, there are nine classes in LAPACK.h++ for solving least squares problems. The approach is to construct a factorization of the matrix A that can be used to solve the problem. The classes {TYPE}LeastSq{METHOD} encapsulate this factorization object, where {TYPE} is a float, double, or DComplex, and {METHOD} encodes the type of factorization used. There are three methods provided:
Cholesky method or method of normal equations, encoded by Ch
Complete orthogonal factorization method, encoded by QR
Singular value method, encoded by SV.
Which method you choose depends on the tradeoff between speed and robustness that you want to make. As you move down the list, the methods decrease in speed, but increase in accuracy and robustness.
A system of equations Ax=b is underdetermined if the columns of A are linearly dependent. In this case, there are infinitely many solutions. Often, the one that is desired is the one with minimum norm. The classes {TYPE}LeastSqQR and {TYPE}LeastSqSV allow you to compute this solution.
Note that the classes {TYPE}LeastSqQR and {TYPE}LeastSqSV are publicly derived from their underlying decomposition types, so that advanced users can manipulate the underlying decomposition directly.
// Find the least squares (if overdetermined) or
// minimum norm (if underdetermined) solution
// to the linear system of equations Ax=b. // Also input a tolerance called epsilon which represents the // accuracy of the entries in A. This makes the procedure // robust in the event that A has nearly linearly dependent // columns. #include <iostream.h> #include <rw/dlssv.h> main() { RWGenMat<double> A; double eps; RWMathVec<double> b; cin >> A >> eps >> b; DoubleLeastSqSV decomp(A,tol); RWMathVec<double> x = decomp.solve(b); cout << "Solution is " << x << endl; }
FloatLeastSqCh(); DoubleLeastSqCh(); DComplexLeastSqCh(); FloatLeastSqQR(); DoubleLeastSqQR(); DComplexLeastSqQR(); FloatLeastSqSV(); DoubleLeastSqSV(); DComplexLeastSqSV();
Default constructor. Builds an empty decomposition. The decomposition is filled in using the factor function.
FloatLeastSqCh(const RWGenMat<float>& A); DoubleLeastSqCh(const RWGenMat<double>& A); DComplexLeastSqCh(const RWGenMat<DComplex>& A); FloatLeastSqQR(const RWGenMat<float>& A, float tol=0); DoubleLeastSqQR(const RWGenMat<double>& A, double tol=0); DComplexLeastSqQR(const RWGenMat<DComplex>& A, double tol=0); FloatLeastSqSV(const RWGenMat<float>& A, float tol=0); DoubleLeastSqSV(const RWGenMat<double>& A, double tol=0); DComplexLeastSqSV(const RWGenMat<DComplex>& A, double tol=0);
Constructs a factorization representing the given matrix. In the QR and SV methods, the supplied tolerance tells when to consider entries along the diagonal of the underlying decomposition as 0. The tolerance should indicate the accuracy to which you know entries in A.
FloatLeastSqQR(const FloatQRDecomp& A, float tol=0); DoubleLeastSqQR(const DoubleQRDecomp& A, double tol=0); DComplexLeastSqQR(const DComplexQRDecomp& A, double tol=0); FloatLeastSqQR(const FloatCODecomp& A); DoubleLeastSqQR(const DoubleCODecomp& A); DComplexLeastSqQR(const DComplexCODecomp& A); FloatLeastSqSV(const FloatSVDecomp& A, float tol=0); DoubleLeastSqSV(const DoubleSVDecomp& A, double tol=0); DComplexLeastSqSV(const DComplexSVDecomp& A, double tol=0);
Constructs a factorization representing the linear system represented by the given decomposition. The parameter tol tells when to consider entries along the diagonal of the given QR or SV decomposition as 0. The tolerance should indicate the accuracy to which you know the entries in the matrix.
unsigned FloatLeastSqCh::cols(); unsigned DoubleLeastSqCh::cols(); unsigned DComplexLeastSqCh::cols(); unsigned FloatLeastSqQR::cols(); unsigned DoubleLeastSqQR::cols(); unsigned DComplexLeastSqQR::cols(); unsigned FloatLeastSqSV::cols(); unsigned DoubleLeastSqSV::cols(); unsigned DComplexLeastSqSV::cols();
Returns the number of columns in the matrix represented by this factorization.
void FloatLeastSqCh::factor(const RWGenMat<float>& A); void DoubleLeastSqCh::factor(const RWGenMat<double>& A); void DComplexLeastSqCh::factor(const RWGenMat<DComplex>& A); void FloatLeastSqQR::factor(const RWGenMat<float>& A, float tol=0); void DoubleLeastSqQR::factor(const RWGenMat<double>& A, double tol=0); void DComplexLeastSqQR::factor(const RWGenMat<DComplex>& A, double tol=0); void FloatLeastSqSV::factor(const RWGenMat<float>& A, float tol=0); void DoubleLeastSqSV::factor(const RWGenMat<double>& A, double tol=0); void DComplexLeastSqSV::factor(const RWGenMat<DComplex>& A, double tol=0);
Replaces the current factorization with a factorization representing the given matrix. The current factorization is lost. In the QR and SV methods, the supplied tolerance tells when to consider entries along the diagonal of the underlying decomposition as 0. The tolerance should indicate the accuracy to which you know entries in A.
void FloatLeastSqQR::factor(const FloatQRDecomp& A, float tol=0); void DoubleLeastSqQR::factor(const DoubleQRDecomp& A,double tol=0); void DComplexLeastSqQR::factor(const DComplexQRDecomp&, double tol=0); void FloatLeastSqQR::factor(const FloatCODecomp& A); void DoubleLeastSqQR::factor(const DoubleCODecomp& A); void DComplexLeastSqQR::factor(const DComplexCODecomp& A); void FloatLeastSqSV::factor(const FloatSVDecomp& A, float tol=0); void DoubleLeastSqSV::factor(const DoubleSVDecomp& A, double tol=0); void DComplexLeastSqSV::factor(const DComplexSVDecomp&, double tol=0);
Replaces the current factorization with a factorization representing the linear system represented by the given decomposition. The current factorization is lost. The parameter tol tells when to consider entries along the diagonal of the given QR or SV decomposition as 0. The tolerance should indicate the accuracy to which you know the entries in the matrix.
unsigned FloatLeastSqCh::fail(); unsigned DoubleLeastSqCh::fail(); unsigned DComplexLeastSqCh::fail(); unsigned FloatLeastSqSV::fail(); unsigned DoubleLeastSqSV::fail(); unsigned DComplexLeastSqSV::fail(); unsigned FloatLeastSqCh::good(); unsigned DoubleLeastSqCh::good(); unsigned DComplexLeastSqCh::good(); unsigned FloatLeastSqSV::good(); unsigned DoubleLeastSqSV::good(); unsigned DComplexLeastSqSV::good();
Indicates whether the construction of the factorization is successful. In the Cholesky case, the factorization can fail if the system matrix is not of full rank. In the SV case, the factorization can fail if some of the singular values fail to converge.
unsigned FloatLeastSqCh::rank(); unsigned DoubleLeastSqCh::rank(); unsigned DComplexLeastSqCh::rank(); unsigned FloatLeastSqQR::rank(); unsigned DoubleLeastSqQR::rank(); unsigned DComplexLeastSqQR::rank(); unsigned FloatLeastSqSV::rank(); unsigned DoubleLeastSqSV::rank(); unsigned DComplexLeastSqSV::rank();
Returns the rank the matrix represented by this factorization. If a tolerance is used to construct the factorization, this is a measure of the numerical rank of the matrix, that is, the rank assuming entries on the diagonal of the decomposition less than the tolerance are 0.
unsigned FloatLeastSqCh::rows(); unsigned DoubleLeastSqCh::rows(); unsigned DComplexLeastSqCh::rows(); unsigned FloatLeastSqQR::rows(); unsigned DoubleLeastSqQR::rows(); unsigned DComplexLeastSqQR::rows(); unsigned FloatLeastSqSV::rows(); unsigned DoubleLeastSqSV::rows(); unsigned DComplexLeastSqSV::rows();
Returns the number of rows in the matrix represented by this factorization.
void FloatLeastSqSV::truncate(float tol); unsigned DoubleLeastSqSV::truncate(double tol); unsigned DComplexLeastSqSV::truncate(double tol);
Truncates the factorization so that singular values less than the indicated tolerance are treated as 0.
unsigned FloatLeastSqCh::residual(const RWMathVec<float>& b); unsigned DoubleLeastSqCh::residual(const RWMathVec<double>& b); unsigned DComplexLeastSqCh::residual(const RWMathVec<DComplex>& b); unsigned FloatLeastSqQR::residual(const RWMathVec<float>& b); unsigned DoubleLeastSqQR::residual(const RWMathVec<double>& b); unsigned DComplexLeastSqQR::residual(const RWMathVec<DComplex>& b); unsigned FloatLeastSqSV::residual(const RWMathVec<float>& b); unsigned DoubleLeastSqSV::residual(const RWMathVec<double>& b); unsigned DComplexLeastSqSV::residual(const RWMathVec<DComplex>& b);
Returns the residual vector, b-Ax, associated with the right-hand side vector b and the corresponding least squares solution vector x.
unsigned FloatLeastSqCh::residualNorm(const RWMathVec<float>& b); unsigned DoubleLeastSqCh::residualNorm(const RWMathVec<double>& b); unsigned DComplexLeastSqCh::residualNorm(const RWMathVec<DComplex>& b); unsigned FloatLeastSqQR::residualNorm(const RWMathVec<float>& b); unsigned DoubleLeastSqQR::residualNorm(const RWMathVec<double>& b); unsigned DComplexLeastSqQR::residualNorm(const RWMathVec<DComplex>& b); unsigned FloatLeastSqSV::residualNorm(const RWMathVec<float>& b); unsigned DoubleLeastSqSV::residualNorm(const RWMathVec<double>& b); unsigned DComplexLeastSqSV::residualNorm(const RWMathVec<DComplex>& b);
Returns the norm of the residual vector, (b-Ax)'(b-Ax), associated with the right-hand side vector b and the corresponding least squares solution vector x.
RWMathVec<float> FloatLeastSqCh::solve(const RWMathVec<float>& b); RWMathVec<double> DoubleLeastSqCh::solve(const RWMathVec<double>& b); RWMathVec<DComplex> DComplexLeastSqCh::solve(const RWMathVec<DComplex>& b); RWMathVec<float> FloatLeastSqQR::solve(const RWMathVec<float>& b); RWMathVec<double> DoubleLeastSqQR::solve(const RWMathVec<double>& b); RWMathVec<DComplex> DComplexLeastSqQR::solve(const RWMathVec<DComplex>& b); RWMathVec<float> FloatLeastSqSV::solve(const RWMathVec<float>& b); RWMathVec<double> DoubleLeastSqSV::solve(const RWMathVec<double>& b); RWMathVec<DComplex> DComplexLeastSqSV::solve(const RWMathVec<DComplex>& b);
Returns the solution to the least squares problem for the given right-hand-side vector b. If the system matrix is overdetermined, this is the vector that minimizes the norm of the residual vector. If the system matrix is not overdetermined, this is the minimum norm vector that satisfies the system equations.
©Copyright 1999, Rogue Wave Software, Inc.
Send mail to report errors or comment on the documentation.