synaps/upol/EUCLIDIAN.h

#ifndef SYNAPS_UPOL_EUCLIDIAN_H
#define SYNAPS_UPOL_EUCLIDIAN_H

#include <synaps/init.h>
#include <synaps/base/type.h>
#include <synaps/arithm/pow.h>
#include <synaps/upol/UPOLDAR.m>


__BEGIN_NAMESPACE_SYNAPS


namespace EUCLIDIAN
{
  using UPOLDAR::lcoeff;
  using UPOLDAR::degree;

    template<class R>
    void div_rem(R & q, R & a, const R & b)
    {

      typedef typename R::iterator iterator;
      R t;
      typename R::value_type lb=b[degree(b)];
        q=0;
        while (degree(a)>=degree(b)) {
          typename R::value_type p = a[degree(a)]/lb;
          R m( a[degree(a)]/lb, degree(a)-degree(b) );
          t=b*m;
          t[degree(a)]=a[degree(a)];
          a-=t;
          q += m;
        }
    }
    
  template<class R>
        void pseudo_div_rem(R & q, R & a, const R & b)
        {
          typedef typename R::iterator iterator;
          R t;
          typename R::value_type lb=b[degree(b)];
          q=0;
          while (degree(a)>=degree(b)) {
                R m(a[degree(a)],degree(a)-degree(b));
                t=b*m;
                a*=lb;
                t[degree(a)]=a[degree(a)];
                a-=t;
                q+=m;
          }
        }


  template<class Poly>
        void prem(const Poly& f, const Poly& g, Poly& q, Poly & r)
        {
          r = f * pow(lcoeff(g), degree(f) - degree(g) + 1);
          q = 0;
          while (degree(r) >= degree(g)) {
                int delta =  degree(r) - degree(g);
                q = q + Poly(1,delta) * (lcoeff(r)/lcoeff(g));
                r = r - ( (g * Poly(1, delta)) * (lcoeff(r)/lcoeff(g)) );
          }
          return;
        }


  template<class Poly>
        Poly prem(const Poly& f, const Poly& g, Poly& q)
        {
          Poly r(f * pow(lcoeff(g), degree(f) - degree(g) + 1) );
          q = 0;
          while (degree(r) >= degree(g)) {
                int delta =  degree(r) - degree(g);
                q = q + Poly(1,delta) * (lcoeff(r)/lcoeff(g));
                r = r - ( (g * Poly(1, delta)) * (lcoeff(r)/lcoeff(g)) );
          }
          return r;
        }


  template<class Poly>
    Poly prem(const Poly& f, const Poly& g)
    {
      Poly r(f * pow(lcoeff(g), degree(f) - degree(g) + 1) );
      Poly q = 0;
      while (degree(r) >= degree(g)) {
        int delta =  degree(r) - degree(g);
        //              q = q + Poly(1,delta) * (lcoeff(r)/lcoeff(g));
        r = r - ( (g * Poly(1, delta)) * (lcoeff(r)/lcoeff(g)) );
      }
      return r;
    }
  
  template<class Poly>
    Poly pquo(const Poly& f, const Poly& g)
    {
      Poly r(f * pow(lcoeff(g), degree(f) - degree(g) + 1) );
      Poly q = 0;
      while (degree(r) >= degree(g)) {
        int delta =  degree(r) - degree(g);
        q = q + Poly(1,delta) * (lcoeff(r)/lcoeff(g));
        r = r - ( (g * Poly(1, delta)) * (lcoeff(r)/lcoeff(g)) );
      }
      return q;
    }
  
  
  
  template<class R>
    void div_pseudorem0(R & q, R & a, const R & b)
    {
      typedef typename R::iterator iterator;
      R t(degree(b)+1, type::AsSize());
      typename R::value_type lb=b[degree(b)];
      q=0;
      int c=0, dab =degree(a)-degree(b)+1;
      while (degree(a)>=degree(b)) {
        R m(a[degree(a)],degree(a)-degree(b));
        t=b*m;
        c++;
        a*=lb;
        t[degree(a)]=a[degree(a)];a-=t;
        q+=m;
      }
      a *= pow(lb, dab-c);
    }
  
}
//CLOSE


__END_NAMESPACE_SYNAPS

#endif // SYNAPS_UPOL_EUCLIDIAN