powerseries_division_algorithm()

Nterm * PolyRing::powerseries_division_algorithm ( Nterm * f, Nterm * g, Nterm *& quot ) const

protected

Definition at line 1553 of file poly.cpp.

{
  // This is intended for use when there is one variable, inverses are present,
  // and the coefficient ring is a field, or is ZZ.
  // The algorithm used is as follows.
  // Given a non-zero polynomial f = f(t,t^-1), let v(f) = top - bottom
  //   where top is the largest exponent in f, and bottom is the smallest.
  //   So if f is a monomial, v(f) = 0.  Also, v(fg) = v(f)+v(g) (at least if
  //   the
  //   coefficient ring is a domain), and v(f) >= 0 always.
  // The algorithm is as follows:
  //   Reduce f = f0 by lt(g) to obtain f1, then again to obtain f2, etc.
  //   So v(f1) >= v(f2) >= ... >= v(fi),
  //   and either fi = 0, v(fi) < v(g), or v(f(i+1)) > v(fi).
  //   In this case, the remainder returned is fi.
  //   (Note: the last case won't happen if the coefficients are a field, or the
  //   lead coefficient of g is a unit).
 
  // This is intended for use when the monomial order has an initial weight
  // vector and the first
  //   term of the divisor's weight value is strictly greater than the other
  //   terms.
  //   AND where (all) inverses are present,
  // and the coefficient ring is a field, or is ZZ.
  // The algorithm used is as follows.
  // Given a non-zero polynomial f = f(t,t^-1), let v(f) = top - bottom
  //   where top is the weight vec value of the first monomial of f, and bottom
  //   is the weight vec value
  //   for the last term of f.
  //   So if f is a monomial, v(f) = 0.  Also, v(fg) = v(f)+v(g) (at least if
  //   the
  //   coefficient ring is a domain), and v(f) >= 0 always.
  // The algorithm is as follows:
  //   Reduce f = f0 by lt(g) to obtain f1, then again to obtain f2, etc.
  //   So v(f1) >= v(f2) >= ... >= v(fi),
  //   and either fi = 0, v(fi) < v(g), or v(f(i+1)) > v(fi).
  //   In this case, the remainder returned is fi.
  //   (Note: the last case won't happen if the coefficients are a field, or the
  //   lead coefficient of g is a unit).
 
  // This returns the remainder, and sets quot to be the quotient.
 
  ring_elem a = copy(f);
  Nterm *t = a;
  Nterm *b = g;
 
  if (!M_->weight_value_exists())
    {
      throw exc::engine_error("no method for division in this ring");
    }
 
  Nterm divhead;
  Nterm remhead;
  Nterm *divt = &divhead;
  Nterm *remt = &remhead;
  Nterm *nextterm = new_term();
  long gval = 0, flast = 0;
 
  if (t != nullptr)
    {
      Nterm *z = t;
      for (; z->next != nullptr; z = z->next)
        ;
 
      flast = M_->first_weight_value(z->monom);
    }
 
  if (b != nullptr)
    {
      long gfirst, glast;
      Nterm *z = b;
 
      if (z->next == nullptr)
        gval = 0;
      else
        {
          gfirst = M_->first_weight_value(z->monom);
          long gsecond = M_->first_weight_value(z->next->monom);
          if (gfirst == gsecond)
            {
              // In this case, we return silently
              throw exc::internal_error("division algorithm for this element is not implemented");
            }
          for (; z->next != nullptr; z = z->next)
            ;
          glast = M_->first_weight_value(z->monom);
          gval = gfirst - glast;
        }
    }
 
  //  buffer o;
  while (t != nullptr)
    {
      long ffirst;
 
      ffirst = M_->first_weight_value(t->monom);
      long fval = ffirst - flast;
 
      // std::cout << "f: "; dNterm(this, t); std::cout << std::endl;
      // std::cout << "g: "; dNterm(this, g); std::cout << std::endl;
      if (fval >= gval and M_->divides(g->monom, t->monom))
        //if (fval >= gval)
        {
          // o << "t = "; elem_text_out(o,t); o << newline;
          a = t;
          bool cancelled = imp_attempt_to_cancel_lead_term(
              a, g, nextterm->coeff, nextterm->monom);
          t = a;
          //    o << "  new t = "; elem_text_out(o,t); o << newline;
          //      o << "  cancelled = " << (cancelled ? "true" : "false") <<
          //      newline;
          //    o << "  coeff = "; K_->elem_text_out(o,nextterm->coeff); o <<
          //    newline;
          //    emit(o.str());
          if (!K_->is_zero(nextterm->coeff))
            {
              divt->next = nextterm;
              divt = divt->next;
              nextterm = new_term();
            }
          if (!cancelled)
            {
              remt->next = t;
              remt = remt->next;
              t = t->next;
            }
        }
      else
        {
          remt->next = t;
          remt = remt->next;
          t = t->next;
        }
    }
 
  nextterm = nullptr;
  remt->next = nullptr;
  divt->next = nullptr;
  quot = divhead.next;
  return remhead.next;
}

References Nterm::coeff, copy(), imp_attempt_to_cancel_lead_term(), PolynomialRing::K_, PolynomialRing::M_, Nterm::monom, new_term(), and Nterm::next.

Referenced by remainderAndQuotient().