polynomial_root.h

Go to the documentation of this file.
00001 /* === S Y N F I G ========================================================= */
00021 /* ========================================================================= */
00022 
00023 /* === S T A R T =========================================================== */
00024 
00025 #ifndef __SYNFIG_POLYNOMIAL_ROOT_H
00026 #define __SYNFIG_POLYNOMIAL_ROOT_H
00027 
00028 /* === H E A D E R S ======================================================= */
00029 
00030 #include <complex>
00031 #include <vector>
00032 
00033 /* === M A C R O S ========================================================= */
00034 
00035 /* === T Y P E D E F S ===================================================== */
00036 
00037 /* === C L A S S E S & S T R U C T S ======================================= */
00038 template < typename T = float, typename F = float >
00039 class Polynomial : public std::vector<T> //a0 + a1x + a2x^2 + ... + anx^n
00040 {
00041 public:
00042 
00043     //Will maintain all lower constants
00044     void degree(unsigned int d, const T & def = (T)0) { resize(d+1,def); }
00045     unsigned int degree()const { return size() - 1; }
00046 
00047     const Polynomial & operator+=(const Polynomial &p)
00048     {
00049         if(p.size() > size())
00050             resize(p.size(), (T)0);
00051 
00052         for(int i = 0; i < p.size(); ++i)
00053         {
00054             (*this)[i] += p[i];
00055         }
00056         return *this;
00057     }
00058 
00059     const Polynomial & operator-=(const Polynomial &p)
00060     {
00061         if(p.size() > size())
00062             resize(p.size(), (T)0);
00063 
00064         for(int i = 0; i < p.size(); ++i)
00065         {
00066             (*this)[i] -= p[i];
00067         }
00068         return *this;
00069     }
00070 
00071     const Polynomial & operator*=(const Polynomial &p)
00072     {
00073         if(p.size() < 1)
00074         {
00075             resize(0);
00076             return *this;
00077         }
00078 
00079         unsigned int i,j;
00080         std::vector<T> nc(*this);
00081 
00082         //in place for constant stuff
00083         for(i = 0; i < nc.size(); ++i)
00084         {
00085             (*this)[i] *= p[0];
00086         }
00087 
00088         if(p.size() < 2) return *this;
00089 
00090         resize(size() + p.degree());
00091         for(int i = 0; i < nc.size(); ++i)
00092         {
00093             for(int j = 1; j < p.size(); ++j)
00094             {
00095                 nc[i+j] += nc[i]*p[j];
00096             }
00097         }
00098 
00099         return *this;
00100     }
00101 };
00102 
00103 class RootFinder
00104 {
00105     std::vector< std::complex<float> >  workcoefs;
00106     int its;
00107 
00108 public:
00109     std::vector< std::complex<float> >  coefs; //the number of coefficients determines the degree of polynomial
00110 
00111     std::vector< std::complex<float> >  roots;
00112 
00113     void find_all_roots(bool polish);
00114 };
00115 
00116 
00117 
00118 /* === E N D =============================================================== */
00119 
00120 #endif

Generated on Wed Aug 15 05:00:23 2007 for synfig by  doxygen 1.5.3