C++ Polynomial::degree方法代码示例

void Polynomial<int>::pseudo_div(
  const Polynomial<int>& f, const Polynomial<int>& g, 
  Polynomial<int>& q, Polynomial<int>& r, int& D)
{
  CGAL_NEF_TRACEN("pseudo_div "<<f<<" , "<< g);
  int fd=f.degree(), gd=g.degree();
  if ( fd<gd ) 
  { q = Polynomial<int>(0); r = f; D = 1; 
    CGAL_postcondition(Polynomial<int>(D)*f==q*g+r); return; 
  }
  // now we know fd >= gd and f>=g
  int qd=fd-gd, delta=qd+1, rd=fd;
  { q = Polynomial<int>( std::size_t(delta) ); }; // workaround for SUNPRO
  int G = g[gd]; // highest order coeff of g
  D = G; while (--delta) D*=G; // D = G^delta
  Polynomial<int> res = Polynomial<int>(D)*f;
  CGAL_NEF_TRACEN("  pseudo_div start "<<res<<" "<<qd<<" "<<q.degree());
  while (qd >= 0) {
    int F = res[rd]; // highest order coeff of res
    int t = F/G;     // ensured to be integer by multiplication of D
    q.coeff(qd) = t;    // store q coeff
    res.minus_offsetmult(g,t,qd); 
    if (res.is_zero()) break;
    rd = res.degree();
    qd = rd - gd;
  }
  r = res;
  CGAL_postcondition(Polynomial<int>(D)*f==q*g+r);
  CGAL_NEF_TRACEN("  returning "<<q<<", "<<r<<", "<< D);
}

Polynomial::Polynomial(const Polynomial& rhs)
{
  mMonomials.resize(rhs.degree() + 1, Monomial("0", 0));
  for (int i = 0; i <= rhs.degree(); i++)
    {
      mMonomials[i] = rhs[i];
    }
}

Polynomial Polynomial::operator *(Polynomial b) {
	Polynomial ret(degree() + b.degree() - 1);
	for (int i = 0; i < degree(); i++) {
		for (int j = 0; j < b.degree(); j++) {
			ret.coefficients[i + j] += (coefficients[i] * b.coefficients[j]);
		}
	}
	return ret;
}

Polynomial Polynomial::operator+(Polynomial b) {
	Polynomial ret = *this;
	for (int i = 0; i < std::max(this->degree(), b.degree()); i++) {
		if (i > this->degree()) {
			ret.coefficients.push_back(b.coefficients[i]);
		}
		ret.coefficients[i] += i < b.degree() ? b.coefficients[i] : 0;
	}
	return ret;
}

Polynomial
Polynomial::operator*(const Polynomial& other) const
{
  Polynomial result;

  for (int i = 0; i <= degree() + other.degree(); i++)
    {
      for (int j = MAX(0, i - degree()); 
	   j <= MIN(other.degree(), i);
	   j++)
	{
	  result += mMonomials[i-j] * other.mMonomials[j];	  
	}
    }  

  return result;
}

inline Polynomial<eT>
gcd(const Polynomial<eT>& x, const Polynomial<eT>& y)
{
    Polynomial<eT> r;
    Polynomial<eT> u;
    Polynomial<eT> v;

    const int xdeg = x.degree();
    const int ydeg = y.degree();

    if(xdeg < 0 || ydeg < 0)
        return std::move(u);

    if((xdeg == 0 && x[0] == 0) || (ydeg == 0 && y[0] == 0))
    {
        u.resize(1);
        u[0] = eT(0);
        return std::move(u);
    }

    const bool maxo = (xdeg >= ydeg);
    u = maxo ? x : y;
    v = maxo ? y : x;

    while(1)
    {
        r = u % v;
        u = v;
        v = r;
        if(v.degree() == 0 && v[0] == 0)
            break;
    }

    if(u.degree() == 0) //x and y are co-primes
        u[0] = eT(1);
    return std::move(u);
}

// doesn't return a remainder
// N.B. this function will break if there is a remainder
const Polynomial Polynomial::div( const Polynomial& rhs ) const {
  Polynomial retVal;

  if( degree() < rhs.degree() )
    {
      return retVal; // return 0
    }

  Polynomial rSide( *this );
  int rDeg = rhs.degree();
  double rCoeff = rhs._list.begin().getData().coeff;
  
  itr_t it = rSide._list.begin();
  while( 1 )
    {
      if( it == rSide._list.end() ) break;
      int deg_diff = it.getData().degree() - rDeg;
      if( deg_diff < 0 ) break; // TODO: check this condition, maybe need to put rest into remainder?
      
      double coeff = it.getData().coeff / rCoeff;
      Polynomial tmp;
      Term multiplier( coeff, deg_diff );
      retVal.addTerm( multiplier );
      
      for( itr_t itt = rhs._list.begin(); itt != rhs._list.end(); ++itt )
	{
	  Term res = itt.getData() * multiplier;
	  tmp.addTerm( res );
	}

      rSide = rSide.sub( tmp );
      it = rSide._list.begin();
    }
  
  return retVal;
}

Polynomial Polynomial::operator*(const Polynomial& other) const
{
	Polynomial result;
	int n = (other.degree()) + degree();
	double c = 0;
	int i = 0;
	
	for(int left = 0; left < size; left++){
		for(int right = 0; right < other.size; right++){
			c = coeff[left] * other.coeff[right];
			i = left + right;
			result[i] += c;
		}
	}
	result.size = n + 1;
	return result;
}

inline void
poly2mat(const Polynomial<eT>& p, Mat<eT>& m)
{
    const uword n = p.degree();
    m.set_size(n,n);

    eT head = p[n];

    for(uword i = 0; i < n-1; ++i)
    {
        m(i, i+1) = 1.0;
    }
    for(uword i = 0; i < n; ++i)
    {
        m(n-1, i) = -p[i]/head;
    }
}

FiniteField::FiniteField(const Polynomial& irred) : p(irred.getMod()), k(irred.degree()), domain(irred.getDomain()), modulus(irred) {
    //assert(k > 1 && "Use PrimeField if k==1");
    assert(irred.isIrreducible());
}