Main Page   Groups   Namespace List   Class Hierarchy   Alphabetical List   Compound List   File List   Namespace Members   Compound Members   File Members   Concepts

vnl_symmetric_eigensystem.h

Go to the documentation of this file.
00001 #ifndef vnl_symmetric_eigensystem_h_
00002 #define vnl_symmetric_eigensystem_h_
00003 
00004 //:
00005 //  \file
00006 //  \brief Find eigenvalues of a symmetric matrix
00007 //  \author Andrew W. Fitzgibbon, Oxford RRG, 29 Aug 96
00008 //  \verbatim
00009 //  Modifications
00010 //  fsm@robots, 5 March 2000: templated
00011 //  dac (Manchester) 28/03/2001: tidied up documentation
00012 //  \endverbatim
00013 //  
00014 
00015 #include <vnl/vnl_matrix.h>
00016 #include <vnl/vnl_diag_matrix.h>
00017 
00018 
00019 //: Find eigenvalues of a symmetric matrix
00020 //
00021 // 
00022 //    Solve the eigenproblem \f$A x = \lambda x\f$, with \f$A\f$ symmetric.
00023 //    The resulting eigenvectors and values are sorted in increasing order
00024 //    so <CODE> V.column(0) </CODE> is the eigenvector corresponding to the smallest
00025 //    the smallest eigenvalue.
00026 //
00027 //    As a matrix decomposition, this is \f$A = V D V^t\f$
00028 //
00029 //    Uses the EISPACK routine RS, which in turn calls TRED2 to reduce A
00030 //    to tridiagonal form, followed by TQL2, to find the eigensystem.
00031 //    This is summarized in Golub and van Loan, \S8.2.  The following are
00032 //    the original subroutine headers:
00033 //
00034 // \remark TRED2 is a translation of the Algol procedure tred2,
00035 //     Num. Math. 11, 181-195(1968) by Martin, Reinsch, and Wilkinson.
00036 //     Handbook for Auto. Comp., Vol.ii-Linear Algebra, 212-226(1971).
00037 //
00038 // \remark This subroutine reduces a real symmetric matrix to a
00039 //     symmetric tridiagonal matrix using and accumulating
00040 //     orthogonal similarity transformations.
00041 //
00042 // \remark TQL2 is a translation of the Algol procedure tql2,
00043 //     Num. Math. 11, 293-306(1968) by Bowdler, Martin, Reinsch, and Wilkinson.
00044 //     Handbook for Auto. Comp., Vol.ii-Linear Algebra, 227-240(1971).
00045 //
00046 // \remark This subroutine finds the eigenvalues and eigenvectors
00047 //     of a symmetric tridiagonal matrix by the QL method.
00048 //     the eigenvectors of a full symmetric matrix can also
00049 //     be found if  tred2  has been used to reduce this
00050 //     full matrix to tridiagonal form.
00051 
00052 bool vnl_symmetric_eigensystem_compute(vnl_matrix<float> const & A, 
00053                                        vnl_matrix<float> & V,
00054                                        vnl_vector<float> & D);
00055 
00056 //: Find eigenvalues of a symmetric matrix
00057 //
00058 // 
00059 //    Solve the eigenproblem \f$A x = \lambda x\f$, with \f$A\f$ symmetric.
00060 //    The resulting eigenvectors and values are sorted in increasing order
00061 //    so <CODE> V.column(0) </CODE> is the eigenvector corresponding to the smallest
00062 //    the smallest eigenvalue.
00063 bool vnl_symmetric_eigensystem_compute(vnl_matrix<double> const & A, 
00064                                        vnl_matrix<double> & V,
00065                                        vnl_vector<double> & D);
00066 
00067 //: Computes and stores the eigensystem decomposition
00068 // of a symmetric matrix.
00069 export template <class T>
00070 class vnl_symmetric_eigensystem {
00071 public:
00072   //: Solve real symmetric eigensystem \f$A x = \lambda x\f$ 
00073   vnl_symmetric_eigensystem(vnl_matrix<T> const & M);
00074 
00075 protected:
00076   // need this here to get inits in correct order, but still keep gentex
00077   // in the right order.
00078   int n_;
00079 
00080 public:
00081   //: Public eigenvectors.  After construction, the columns of V are the
00082   // eigenvectors, sorted by increasing eigenvalue, from most negative to
00083   // most positive.
00084   vnl_matrix<T> V;
00085 
00086   //: Public eigenvalues.  After construction,  D contains the
00087   // eigenvalues, sorted as described above.  Note that D is a vnl_diag_matrix,
00088   // and is therefore stored as a vcl_vector while behaving as a matrix.
00089   vnl_diag_matrix<T> D;
00090 
00091   //: Recover specified eigenvector after computation.
00092   vnl_vector<T> get_eigenvector(int i) const;
00093 
00094   //: Recover specified eigenvalue after computation.
00095   T             get_eigenvalue(int i) const;
00096 
00097   //: Convenience method to get least-squares nullvector.
00098   // It is deliberate that the signature is the same as on vnl_svd<T>.
00099   vnl_vector<T> nullvector() const { return get_eigenvector(0); }
00100 
00101   //: Return the matrix \f$V  D  V^\top\f$.  This can be useful if you've
00102   // modified \f$D\f$.  So an inverse is obtained using
00103   // \verbatim
00104   //   vnl_symmetric_eigensystem} eig(A);
00105   //   eig.D.invert_in_place}();
00106   //   vnl_matrix<double> Ainverse = eig.recompose();
00107   // \endverbatim
00108   
00109   vnl_matrix<T> recompose() const { return V * D * V.transpose(); }
00110 
00111   //: return the pseudoinverse.
00112   vnl_matrix<T> pinverse() const;
00113 
00114   //: return the square root, if positive semi-definite.
00115   vnl_matrix<T> square_root() const;
00116 
00117   //: return the inverse of the square root, if positive semi-definite.
00118   vnl_matrix<T> inverse_square_root() const;
00119 
00120   //: Solve LS problem M x = b
00121   vnl_vector<T> solve(vnl_vector<T> const & b);
00122 
00123   //: Solve LS problem M x = b
00124   void solve(vnl_vector<T> const & b, vnl_vector<T> * x);
00125 };
00126 
00127 #endif // vnl_symmetric_eigensystem_h_

Generated at Wed Mar 12 01:13:16 2003 for ITK by doxygen 1.2.15 written by Dimitri van Heesch, © 1997-2000