synfig-core/trunk/src/synfig/polynomial_root.h

   1 /* === S Y N F I G ========================================================= */
   2 /*!     \file polynomial_root.h
   3 **      \brief Polynomial Root Finder Header
   4 **
   5 **      $Id: polynomial_root.h,v 1.1.1.1 2005/01/04 01:23:14 darco Exp $
   6 **
   7 **      \legal
   8 **      Copyright (c) 2002 Robert B. Quattlebaum Jr.
   9 **
  10 **      This software and associated documentation
  11 **      are CONFIDENTIAL and PROPRIETARY property of
  12 **      the above-mentioned copyright holder.
  13 **
  14 **      You may not copy, print, publish, or in any
  15 **      other way distribute this software without
  16 **      a prior written agreement with
  17 **      the copyright holder.
  18 **      \endlegal
  19 */
  20 /* ========================================================================= */
  21
  22 /* === S T A R T =========================================================== */
  23
  24 #ifndef __SYNFIG_POLYNOMIAL_ROOT_H
  25 #define __SYNFIG_POLYNOMIAL_ROOT_H
  26
  27 /* === H E A D E R S ======================================================= */
  28
  29 #include <complex>
  30 #include <vector>
  31
  32 /* === M A C R O S ========================================================= */
  33
  34 /* === T Y P E D E F S ===================================================== */
  35
  36 /* === C L A S S E S & S T R U C T S ======================================= */
  37 template < typename T = float, typename F = float >
  38 class Polynomial : public std::vector<T> //a0 + a1x + a2x^2 + ... + anx^n
  39 {
  40 public:
  41
  42         //Will maintain all lower constants
  43         void degree(unsigned int d, const T & def = (T)0) { resize(d+1,def); }
  44         unsigned int degree()const { return size() - 1; }
  45
  46         const Polynomial & operator+=(const Polynomial &p)
  47         {
  48                 if(p.size() > size())
  49                         resize(p.size(), (T)0);
  50
  51                 for(int i = 0; i < p.size(); ++i)
  52                 {
  53                         (*this)[i] += p[i];
  54                 }
  55                 return *this;
  56         }
  57
  58         const Polynomial & operator-=(const Polynomial &p)
  59         {
  60                 if(p.size() > size())
  61                         resize(p.size(), (T)0);
  62
  63                 for(int i = 0; i < p.size(); ++i)
  64                 {
  65                         (*this)[i] -= p[i];
  66                 }
  67                 return *this;
  68         }
  69
  70         const Polynomial & operator*=(const Polynomial &p)
  71         {
  72                 if(p.size() < 1)
  73                 {
  74                         resize(0);
  75                         return *this;
  76                 }
  77
  78                 unsigned int i,j;
  79                 std::vector<T> nc(*this);
  80
  81                 //in place for constant stuff
  82                 for(i = 0; i < nc.size(); ++i)
  83                 {
  84                         (*this)[i] *= p[0];
  85                 }
  86
  87                 if(p.size() < 2) return *this;
  88
  89                 resize(size() + p.degree());
  90                 for(int i = 0; i < nc.size(); ++i)
  91                 {
  92                         for(int j = 1; j < p.size(); ++j)
  93                         {
  94                                 nc[i+j] += nc[i]*p[j];
  95                         }
  96                 }
  97
  98                 return *this;
  99         }
 100 };
 101
 102 class RootFinder
 103 {
 104         std::vector< std::complex<float> >      workcoefs;
 105         int     its;
 106
 107 public:
 108         std::vector< std::complex<float> >      coefs; //the number of coefficients determines the degree of polynomial
 109
 110         std::vector< std::complex<float> >      roots;
 111
 112         void find_all_roots(bool polish);
 113 };
 114
 115
 116
 117 /* === E N D =============================================================== */
 118
 119 #endif