1 #ifndef GINAC_UPOLY_HEUR_GCD
2 #define GINAC_UPOLY_HEUR_GCD
4 #include "ring_traits.hpp"
5 #include "normalize.tcc"
6 #include "remainder.tcc"
8 #include "interpolate_padic_uvar.h"
14 /// Compute GCD of primitive univariate polynomials.
15 template<typename T> static bool
16 heur_gcd_z_pp(T& g, const T& a, const T& b, unsigned max_tries = 66)
18 typedef typename T::value_type ring_t;
19 const ring_t n73794 = get_ring_elt(b[0], 73794);
20 const ring_t n27011 = get_ring_elt(b[0], 27011);
21 const std::size_t maxdeg = std::max(degree(a), degree(b));
26 // find the evaluation point
27 ring_t xi = (std::min(max_coeff(a), max_coeff(b)) + 1) << 1;
30 const ring_t av = eval(a, xi);
31 const ring_t bv = eval(b, xi);
32 const ring_t gamma = gcd(av, bv);
33 interpolate(gg, gamma, xi, maxdeg);
34 normalize_in_ring(gg);
35 remainder_in_ring(r, a, gg);
40 // next evaluation point
41 xi = truncate1(xi*isqrt(isqrt(xi))*n73794, n27011);
42 } while (--max_tries != 0);
46 template<typename T> static bool
47 heur_gcd_z_priv(T& g, const T& a, const T& b, const unsigned max_tries = 66)
49 typedef typename T::value_type ring_t;
52 normalize_in_ring(a_, &acont);
53 normalize_in_ring(b_, &bcont);
54 const ring_t gc = gcd(acont, bcont);
55 bool found = heur_gcd_z_pp(g, a_, b_, max_tries);
65 #endif // GINAC_UPOLY_HEUR_GCD