+ if (a.info(info_flags::integer_polynomial) &&
+ b.info(info_flags::integer_polynomial)) {
+ try {
+ return heur_gcd_z(res, a, b, ca, cb, var);
+ } catch (gcdheu_failed) {
+ return false;
+ }
+ }
+
+ // convert polynomials to Z[X]
+ const numeric a_lcm = lcm_of_coefficients_denominators(a);
+ const numeric ab_lcm = lcmcoeff(b, a_lcm);
+
+ const ex ai = a*ab_lcm;
+ const ex bi = b*ab_lcm;
+ if (!ai.info(info_flags::integer_polynomial))
+ throw std::logic_error("heur_gcd: not an integer polynomial [1]");
+
+ if (!bi.info(info_flags::integer_polynomial))
+ throw std::logic_error("heur_gcd: not an integer polynomial [2]");
+
+ bool found = false;
+ try {
+ found = heur_gcd_z(res, ai, bi, ca, cb, var);
+ } catch (gcdheu_failed) {
+ return false;
+ }
+
+ // GCD is not unique, it's defined up to a unit (i.e. invertible
+ // element). If the coefficient ring is a field, every its element is
+ // invertible, so one can multiply the polynomial GCD with any element
+ // of the coefficient field. We use this ambiguity to make cofactors
+ // integer polynomials.
+ if (found)
+ res /= ab_lcm;
+ return found;
+}
+
+
+// gcd helper to handle partially factored polynomials (to avoid expanding
+// large expressions). At least one of the arguments should be a power.
+static ex gcd_pf_pow(const ex& a, const ex& b, ex* ca, ex* cb, bool check_args);
+
+// gcd helper to handle partially factored polynomials (to avoid expanding
+// large expressions). At least one of the arguments should be a product.
+static ex gcd_pf_mul(const ex& a, const ex& b, ex* ca, ex* cb, bool check_args);
+
+/** Compute GCD (Greatest Common Divisor) of multivariate polynomials a(X)
+ * and b(X) in Z[X]. Optionally also compute the cofactors of a and b,
+ * defined by a = ca * gcd(a, b) and b = cb * gcd(a, b).
+ *
+ * @param a first multivariate polynomial
+ * @param b second multivariate polynomial
+ * @param ca pointer to expression that will receive the cofactor of a, or NULL
+ * @param cb pointer to expression that will receive the cofactor of b, or NULL
+ * @param check_args check whether a and b are polynomials with rational
+ * coefficients (defaults to "true")
+ * @return the GCD as a new expression */
+ex gcd(const ex &a, const ex &b, ex *ca, ex *cb, bool check_args, unsigned options)
+{
+#if STATISTICS
+ gcd_called++;
+#endif