/Users/mourrain/Devel/mmx/algebramix/include/algebramix/algebraic.hpp

Go to the documentation of this file.

/******************************************************************************
* MODULE     : algebraic.hpp
* DESCRIPTION: Algebraic extensions based on a primitive element
* COPYRIGHT  : (C) 2007  Joris van der Hoeven
*******************************************************************************
* This software falls under the GNU general public license and comes WITHOUT
* ANY WARRANTY WHATSOEVER. See the file $TEXMACS_PATH/LICENSE for more details.
* If you don't have this file, write to the Free Software Foundation, Inc.,
* 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA.
******************************************************************************/

#ifndef __MMX_ALGEBRAIC_HPP
#define __MMX_ALGEBRAIC_HPP
#include <algebramix/algebraic_extension.hpp>
namespace mmx {
#define TMPL_DEF template<typename C, typename Extension=algebraic_extension<C> >
#define TMPL template<typename C, typename Extension>
#define Polynomial polynomial<C>
#define Algebraic algebraic<C,Extension>
#define Element typename Extension::El
TMPL class algebraic;
TMPL class primitive_element;

/******************************************************************************
* Trivial extensions
******************************************************************************/

template<typename FT, typename C, typename Extension>
struct trivial_extension_helper {
  static inline Extension ext () {
    return Extension (); }
  static inline Extension ext (const format<C>& fm) {
    return Extension (fm); }
};

template<typename C, typename Extension>
struct trivial_extension_helper<empty_format,C,Extension> {
  static inline Extension ext () {
    static Extension ext= Extension (format<C> ());
    return ext; }
  static inline Extension ext (const format<C>& fm) {
    static Extension ext= Extension (format<C> ());
    return ext; }
};

TMPL inline Extension
trivial_extension () {
  typedef typename format<C>::FT FT;
  return trivial_extension_helper<FT,C,Extension>::ext ();
}

TMPL inline Extension
trivial_extension (const format<C>& fm) {
  typedef typename format<C>::FT FT;
  return trivial_extension_helper<FT,C,Extension>::ext (fm);
}

/******************************************************************************
* The algebraic type
******************************************************************************/

TMPL_DEF
class algebraic {
MMX_ALLOCATORS;
public:
  Extension ext;
  Element p;
public:
  inline algebraic ():
    ext (trivial_extension<C,Extension> ()) {}
  template<typename T> inline algebraic (const T& c):
    ext (trivial_extension<C,Extension> (get_format (as<C> (c)))), p (c) {}
  inline algebraic (const Extension ext2, const Polynomial& p2):
    ext (ext2), p (p2) {}
  inline algebraic (const Polynomial& mp):
    ext (mp), p (vec<C> (promote (0, CF(mp)), promote (1, CF(mp)))) {}
  template<typename X> inline
  algebraic (const Polynomial& mp, const X& x):
    ext (mp, x), p (vec<C> (promote (0, CF(mp)), promote (1, CF(mp)))) {}
  inline algebraic (const Polynomial& mp, const Polynomial& p2):
    ext (mp), p (deg (p2) < deg (mp)? p2: rem (p2, mp)) {}
  template<typename X> inline
  algebraic (const Polynomial& mp, const Polynomial& p2, const X& x):
    ext (mp, x), p (deg (p2) < deg (mp)? p2: rem (p2, mp)) {}
  template<typename C2, typename Extension2> inline
  algebraic (const algebraic<C2,Extension2>& x2):
    ext (field (x2)), p (as<Element > (value (x2))) {}
};

STYPE_TO_TYPE(TMPL,scalar_type,Algebraic,C);

TMPL inline format<C> CF (const Algebraic& x) { return CF(x.p); }
TMPL inline Extension field (const Algebraic& x) { return x.ext; }
TMPL inline Element value (const Algebraic& x) { return x.p; }
TMPL inline Polynomial field_modulus (const Algebraic& x) { return x.ext.mp; }

/******************************************************************************
* Basic routines
******************************************************************************/

TMPL inline nat hash (const Algebraic& x) {
  nat h= hash (value (x));
  return (h<<1) ^ (h<<5) ^ (h>>29) ^ hash (field (x)); }
TMPL inline nat exact_hash (const Algebraic& x) {
  nat h= exact_hash (value (x));
  return (h<<1) ^ (h<<5) ^ (h>>29) ^ exact_hash (field (x)); }
TMPL inline nat hard_hash (const Algebraic& x) {
  nat h= hard_hash (value (x));
  return (h<<1) ^ (h<<5) ^ (h>>29) ^ hard_hash (field (x)); }
TMPL inline bool exact_eq (const Algebraic& x1, const Algebraic& x2) {
  return exact_eq (value (x1), value (x2)) &&
         exact_eq (field (x1), field (x2)); }
TMPL inline bool exact_neq (const Algebraic& x1, const Algebraic& x2) {
  return !exact_eq (x1, x2); }
TMPL inline bool hard_eq (const Algebraic& x1, const Algebraic& x2) {
  return hard_eq (value (x1), value (x2)) &&
         hard_eq (field (x1), field (x2)); }
TMPL inline bool hard_neq (const Algebraic& x1, const Algebraic& x2) {
  return !hard_eq (x1, x2); }

TMPL inline bool is_zero (const Algebraic& x) {
  return is_zero (field (x), value (x)); }
TMPL inline int sign (const Algebraic& x) {
  return sign (field (x), value (x)); }
TMPL inline bool operator == (const Algebraic& x1, const Algebraic& x2) {
  return is_zero (x1 - x2); }
TMPL inline bool operator != (const Algebraic& x1, const Algebraic& x2) {
  return !is_zero (x1 - x2); }

EQUAL_INT_SUGAR(TMPL,Algebraic);
COMPARE_SUGAR(TMPL,Algebraic);
COMPARE_INT_SUGAR(TMPL,Algebraic);

TMPL syntactic
flatten (const Algebraic& x) {
  return apply ("algebraic", flatten (value (x)), flatten (field (x)));
}

TMPL
struct binary_helper<Algebraic >: public void_binary_helper<Algebraic > {
  static inline string short_type_name () {
    return "A" * Short_type_name (C); }
  static inline generic full_type_name () {
    return gen ("Algebraic", Full_type_name (C)); }
  static inline generic disassemble (const Algebraic& x) {
    return gen_vec (as<generic> (field (x)),
                    as<generic> (value (x))); }
  static inline Algebraic assemble (const generic& v) {
    return Algebraic (as<Extension > (vector_access (v, 0)),
                      as<Element > (vector_access (v, 1))); }
  static inline void write (const port& out, const Algebraic& p) {
    binary_write<Extension > (out, field (p));
    binary_write<Element > (out, value (p)); }
  static inline Algebraic read (const port& in) {
    Extension ext= binary_read<Extension > (in);
    Element p= binary_read<Element > (in);
    return Algebraic (ext, p); }
};

/******************************************************************************
* Joining algebraic extensions and minimal annihilators
******************************************************************************/

TMPL inline void
upgrade (Algebraic& a1, Algebraic& a2) {
  if (hard_neq (field (a1), field (a2))) {
    Element p1= value (a1), p2= value (a2);
    Extension ext= upgrade (field (a1), field (a2), p1, p2);
    a1= Algebraic (ext, p1);
    a2= Algebraic (ext, p2);
  }
}

TMPL inline Polynomial
annihilator (const Algebraic& a) {
  return annihilator (field (a), value (a));
}

TMPL inline Algebraic
normalize (const Algebraic& a) {
  Extension  ext= normalize (field (a), value (a));
  Polynomial pri (vec<C> (promote (0, CF(a)), promote (1, CF(a))));
  return Algebraic (ext, pri);
}

/******************************************************************************
* Basic arithmetic
******************************************************************************/

TMPL Algebraic
operator + (const Algebraic& x1b, const Algebraic& x2b) {
  Algebraic x1= x1b, x2= x2b;
  upgrade (x1, x2);
  return Algebraic (field (x1), value (x1) + value (x2));
}

TMPL inline Algebraic
operator + (const Algebraic& x1, const C& x2) {
  return Algebraic (field (x1), value (x1) + x2);
}

TMPL inline Algebraic
operator + (const C& x1, const Algebraic& x2) {
  return Algebraic (field (x2), x1 + value (x2));
}

TMPL inline Algebraic
operator - (const Algebraic& x) {
  return Algebraic (field (x), -value (x));
}

TMPL Algebraic
operator - (const Algebraic& x1b, const Algebraic& x2b) {
  Algebraic x1= x1b, x2= x2b;
  upgrade (x1, x2);
  return Algebraic (field (x1), value (x1) - value (x2));
}

TMPL inline Algebraic
operator - (const Algebraic& x1, const C& x2) {
  return Algebraic (field (x1), value (x1) - x2);
}

TMPL inline Algebraic
operator - (const C& x1, const Algebraic& x2) {
  return Algebraic (field (x2), x1 - value (x2));
}

TMPL inline Algebraic
square (const Algebraic& x1) {
  return Algebraic (field (x1), square (field (x1), value (x1)));
}

TMPL Algebraic
operator * (const Algebraic& x1b, const Algebraic& x2b) {
  Algebraic x1= x1b, x2= x2b;
  upgrade (x1, x2);
  return Algebraic (field (x1), mul (field (x1), value (x1), value (x2)));
}

TMPL inline Algebraic
operator * (const Algebraic& x1, const C& x2) {
  return Algebraic (field (x1), value (x1) * x2);
}

TMPL inline Algebraic
operator * (const C& x1, const Algebraic& x2) {
  return Algebraic (field (x2), x1 * value (x2));
}

TMPL inline Algebraic
invert (const Algebraic& x1) {
  return Algebraic (field (x1), invert (field (x1), value (x1)));
}

TMPL Algebraic
operator / (const C& x1, const Algebraic& x2) {
  return Algebraic (field (x2), div (field (x2), x1, value (x2)));
}

TMPL inline Algebraic
operator / (const Algebraic& x1, const C& x2) {
  return Algebraic (field (x1), value (x1) / x2);
}

TMPL inline Algebraic
operator / (const Algebraic& x1b, const Algebraic& x2b) {
  Algebraic x1= x1b, x2= x2b;
  upgrade (x1, x2);
  return Algebraic (field (x1), div (field (x1), value (x1), value (x2)));
}

ARITH_SCALAR_INT_SUGAR(TMPL,Algebraic);

/******************************************************************************
* Roots
******************************************************************************/

TMPL inline Algebraic
root (const Algebraic& x, const int& p) {
  if (p < 0) return 1 / root (x, -p);
  ASSERT (p != 0, "cannot take zero-th roots");
  Algebraic  y= normalize (x);
  Polynomial z (vec<C> (promote (0, CF(x)), promote (1, CF(x))));
  return Algebraic (ramify (field (y), p), z);
}

TMPL inline Algebraic
sqrt (const Algebraic& x) {
  return root (x, 2);
}

#undef TMPL_DEF
#undef TMPL
#undef Polynomial
#undef Algebraic
#undef Element
} // namespace mmx
#endif // __MMX_ALGEBRAIC_HPP