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

vnl_rpoly_roots.h

Go to the documentation of this file.
00001 #ifndef vnl_rpoly_roots_h_
00002 #define vnl_rpoly_roots_h_
00003 
00004 //:
00005 //  \file
00006 //  \brief Finds roots of a real polynomial
00007 //  \author  Andrew W. Fitzgibbon, Oxford RRG, 06 Aug 96
00008 //
00009 //  Modifications
00010 //  23 may 97, Peter Vanroose - "NO_COMPLEX" option added 
00011 //  (until "complex" type is standardised)
00012 //  dac (Manchester) 28/03/2001: tidied up documentation
00013 
00014 #include <vcl_complex.h>
00015 #include <vnl/vnl_vector.h>
00016 
00017 //: Find the roots of a real polynomial.  
00018 //  Uses algorithm 493 from
00019 //  ACM Trans. Math. Software - the Jenkins-Traub algorithm, described
00020 //  by Numerical Recipes under "Other sure-fire techniques" as
00021 //  "practically a standard in black-box polynomial rootfinders".
00022 //  (See M.A. Jenkins, ACM TOMS 1 (1975) pp. 178-189.).
00023 //
00024 //  This class is not very const-correct as it is intended as a compute object
00025 //  rather than a data object.
00026 
00027 
00028 
00029 class vnl_real_polynomial;
00030 
00031 //: Roots of real polynomial.
00032 class vnl_rpoly_roots {
00033 public:
00034 // Constructors/Destructors--------------------------------------------------
00035 
00036 //: The constructor calculates the roots.  This is the most efficient interface
00037 // as all the result variables are initialized to the correct size.
00038 // @{The polynomial is $ a[0] x^d + a[1] x^{d-1} + \cdots + a[d] = 0 $.@}
00039   vnl_rpoly_roots(const vnl_vector<double>& a);
00040 
00041 //: Calculate roots of RealPolynomial.
00042   vnl_rpoly_roots(const vnl_real_polynomial& poly);
00043 
00044   // Operations----------------------------------------------------------------
00045 
00046 //: Return i'th complex root
00047   vcl_complex<double> operator [] (int i) 
00048         const { return vcl_complex<double>(r_[i], i_[i]); }
00049 
00050 //: Complex vector of all roots.
00051   vnl_vector<vcl_complex<double> > roots() const;
00052 
00053 //: Real part of root I.
00054   const double& real(int i) const { return r_[i]; }
00055 
00056 //: Imaginary part of root I.
00057   const double& imag(int i) const { return i_[i]; }
00058 
00059 //: Vector of real parts of roots
00060   vnl_vector<double>& real() { return r_; }
00061 
00062 //: Vector of imaginary parts of roots
00063   vnl_vector<double>& imag() { return i_; }
00064   
00065 //: Return real roots only.  Roots are real if the absolute value
00066 // of their imaginary part is less than the optional argument TOL.
00067 // TOL defaults to 1e-12 [untested]
00068   vnl_vector<double> realroots(double tol = 1e-12) const;
00069 
00070   // Computations--------------------------------------------------------------
00071 
00072   //: Compute roots using Jenkins-Traub algorithm.
00073   bool compute();
00074   
00075   //: Compute roots using QR decomposition of companion matrix. [unimplemented]
00076   bool compute_qr();
00077   
00078   //: Compute roots using Laguerre algorithm. [unimplemented]
00079   bool compute_laguerre();
00080   
00081   // Data Access---------------------------------------------------------------
00082 
00083   // Data Control--------------------------------------------------------------
00084 
00085 protected:
00086   // Data Members--------------------------------------------------------------
00087   vnl_vector<double> coeffs_;
00088 
00089   vnl_vector<double> r_;
00090   vnl_vector<double> i_;
00091   
00092   int num_roots_found_;
00093   
00094 private:
00095   // Helpers-------------------------------------------------------------------
00096 };
00097 
00098 #endif // vnl_rpoly_roots_h_

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