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

	polynomial operator -(const polynomial& p1, const polynomial& p2)
	{
		size_t j=p1.degree( );
		size_t q= p2.degree( );
		polynomial answer;

		if(j>=q)
		{
			answer.reserve(j+1);
			size_t i;
			for (i=0; i<=j; i++)
			{
				answer.add_to_coef((p1.coefficient(i)-p2.coefficient(i)),(i));
			}
		}
		else
		{
			answer.reserve(q+1);
			size_t i;
			for (i=0; i<=q; i++)
			{
				answer.add_to_coef((p1.coefficient(i)-p2.coefficient(i)),(i));
			}
		}
		return answer;
	}

	static polynomial div_monomial(const polynomial& divident,
			const polynomial& divisor) {
		int_type degree = divident.degree() - divisor.degree();
		int_type coefficient = divident.data.begin()->second
				/ divisor.data.begin()->second;
		return { {degree, coefficient}};
	}

    const polynomial operator+ (polynomial operand) const {
        std::vector<double> result (std::max(degree(), operand.degree()) + 1);

        for (int i = 0; i < result.size(); i++) {
            double x = i <= degree() ? coefficients[i] : 0.0;
            double y = i <= operand.degree() ? operand[i] : 0.0;
            result[i] = x + y;
        }

        return polynomial { result };
    }

    std::pair<polynomial, polynomial> quo_rem (polynomial divisor) const {
        /* Degrees of the remainder and quotient. */
        int kd = divisor.degree();
        int kr = kd - 1;
        int kq = degree() - kd;

        std::vector<double> remainder (kr + 1);
        std::vector<double> quotient (kq + 1);

        auto dbegin = divisor.coefficients.crbegin() + 1;
        auto dend = divisor.coefficients.crend();
        auto qbegin = quotient.end();
        auto qend = quotient.end();

        for (int i = kq; 0 <= i; i--) {
            quotient[i] = (coefficients[kd+i] - ip(dbegin, dend, qbegin, qend)) / divisor[kd];
            --qbegin;
        }

        for (int i = kr; 0 <= i; i--) {
            remainder[i] = coefficients[i] - ip(dbegin, dend, qbegin, qend);
            ++dbegin;
        }

        return { polynomial { quotient }, polynomial { remainder } };
    }

        int descartes_rule(const polynomial<data__> &f, bool positive)
        {
          // catch special case
          if (f.degree()==0)
            return 0;

          // get the coefficients from the polynomial
          std::vector<data__> a(f.degree()+1);
          for (size_t i=0; i<a.size(); ++i)
          {
            a[i]=f.coefficient(i);

            // change the sign of odd powers if want the negative root count
            if (!positive && i%2==1)
              a[i]*=-1;
          }

          return eli::mutil::poly::root::sign_changes(a.begin(), a.end());
        }

std::vector<polynomial> sturm_chain (polynomial p) {
    std::vector<polynomial> chain (p.degree() + 1);

    chain[0] = p;
    chain[1] = p.derivative();

    for (int i = 2; i < chain.size(); i++) {
        chain[i] = -(chain[i-2] % chain[i-1]);
    }

    return chain;
}

        int sturm_count(const polynomial<data__> &f, const data2__ &xmin, const data2__ &xmax)
        {
          // short circuit degenerate cases
          if (xmax<=xmin)
            return 0;

          // catch another degenerate case
          if (f.degree()==0)
            return 0;

          std::vector< polynomial<data__> > sturm_fun;

          // calculate the Sturm functions
          sturm_functions(sturm_fun, f);

          return sturm_count(sturm_fun, xmin, xmax);
        }

        void sturm_functions(std::vector< polynomial<data__> > &sturm_fun, const polynomial<data__> &f)
        {
          size_t i, n;

          // resize the Sturm function vector
          n=f.degree()+1;
          sturm_fun.resize(n);

          polynomial<data__> *ptemp(f.fp()), pjunk;

          // initialize the Sturm functions
          sturm_fun[0]=f;
          ptemp=f.fp();
          sturm_fun[1]=*ptemp;
          delete ptemp;

          // build the Sturm functions
          sturm_fun[0].divide(std::abs(sturm_fun[0].coefficient(sturm_fun[0].degree())));
          sturm_fun[1].divide(std::abs(sturm_fun[1].coefficient(sturm_fun[1].degree())));
          for (i=2; i<n; ++i)
          {
            if (sturm_fun[i-1].degree()>0)
            {
              data__ small_no, max_c(0);
              typename polynomial<data__>::coefficient_type c;

              sturm_fun[i-1].get_coefficients(c);
              for (typename polynomial<data__>::index_type j=0; j<=sturm_fun[i-1].degree(); ++j)
              {
                if (std::abs(c(j))>max_c)
                  max_c=std::abs(c(j));
              }
              small_no=max_c*std::sqrt(std::numeric_limits<data__>::epsilon());

              pjunk.divide(sturm_fun[i], sturm_fun[i-2], sturm_fun[i-1]);
              sturm_fun[i].negative();
              if (std::abs(sturm_fun[i].coefficient(sturm_fun[i].degree()))>small_no)
                sturm_fun[i].divide(std::abs(sturm_fun[i].coefficient(sturm_fun[i].degree())));
              sturm_fun[i].adjust_zero(small_no);
            }
            else
            {
              sturm_fun[i]=0;
            }
          }
        }

void view(const polynomial& test)
{
    cout << test
	 << " (degree is " << test.degree( ) << ")" << endl;
}