ex dummy;
if (divide_in_z(p, g, ca ? *ca : dummy, var) && divide_in_z(q, g, cb ? *cb : dummy, var)) {
g *= gc;
- ex lc = g.lcoeff(x);
- if (is_exactly_a<numeric>(lc) && ex_to<numeric>(lc).is_negative())
- return -g;
- else
- return g;
+ return g;
}
}
// Input polynomials of the form poly^n are sometimes also trivial
if (is_exactly_a<power>(a)) {
ex p = a.op(0);
+ const ex& exp_a = a.op(1);
if (is_exactly_a<power>(b)) {
- if (p.is_equal(b.op(0))) {
+ ex pb = b.op(0);
+ const ex& exp_b = b.op(1);
+ if (p.is_equal(pb)) {
// a = p^n, b = p^m, gcd = p^min(n, m)
- ex exp_a = a.op(1), exp_b = b.op(1);
if (exp_a < exp_b) {
if (ca)
*ca = _ex1;
*cb = _ex1;
return power(p, exp_b);
}
- }
+ } else {
+ ex p_co, pb_co;
+ ex p_gcd = gcd(p, pb, &p_co, &pb_co, check_args);
+ if (p_gcd.is_equal(_ex1)) {
+ // a(x) = p(x)^n, b(x) = p_b(x)^m, gcd (p, p_b) = 1 ==>
+ // gcd(a,b) = 1
+ if (ca)
+ *ca = a;
+ if (cb)
+ *cb = b;
+ return _ex1;
+ // XXX: do I need to check for p_gcd = -1?
+ } else {
+ // there are common factors:
+ // a(x) = g(x)^n A(x)^n, b(x) = g(x)^m B(x)^m ==>
+ // gcd(a, b) = g(x)^n gcd(A(x)^n, g(x)^(n-m) B(x)^m
+ if (exp_a < exp_b) {
+ return power(p_gcd, exp_a)*
+ gcd(power(p_co, exp_a), power(p_gcd, exp_b-exp_a)*power(pb_co, exp_b), ca, cb, false);
+ } else {
+ return power(p_gcd, exp_b)*
+ gcd(power(p_gcd, exp_a - exp_b)*power(p_co, exp_a), power(pb_co, exp_b), ca, cb, false);
+ }
+ } // p_gcd.is_equal(_ex1)
+ } // p.is_equal(pb)
+
} else {
if (p.is_equal(b)) {
// a = p^n, b = p, gcd = p
if (cb)
*cb = _ex1;
return p;
+ }
+
+ ex p_co, bpart_co;
+ ex p_gcd = gcd(p, b, &p_co, &bpart_co, false);
+
+ if (p_gcd.is_equal(_ex1)) {
+ // a(x) = p(x)^n, gcd(p, b) = 1 ==> gcd(a, b) = 1
+ if (ca)
+ *ca = a;
+ if (cb)
+ *cb = b;
+ return _ex1;
+ } else {
+ // a(x) = g(x)^n A(x)^n, b(x) = g(x) B(x) ==> gcd(a, b) = g(x) gcd(g(x)^(n-1) A(x)^n, B(x))
+ return p_gcd*gcd(power(p_gcd, exp_a-1)*power(p_co, exp_a), bpart_co, ca, cb, false);
}
- }
+ } // is_exactly_a<power>(b)
+
} else if (is_exactly_a<power>(b)) {
ex p = b.op(0);
if (p.is_equal(a)) {
*cb = power(p, b.op(1) - 1);
return p;
}
+
+ ex p_co, apart_co;
+ const ex& exp_b(b.op(1));
+ ex p_gcd = gcd(a, p, &apart_co, &p_co, false);
+ if (p_gcd.is_equal(_ex1)) {
+ // b=p(x)^n, gcd(a, p) = 1 ==> gcd(a, b) == 1
+ if (ca)
+ *ca = a;
+ if (cb)
+ *cb = b;
+ return _ex1;
+ } else {
+ // there are common factors:
+ // a(x) = g(x) A(x), b(x) = g(x)^n B(x)^n ==> gcd = g(x) gcd(g(x)^(n-1) A(x)^n, B(x))
+
+ return p_gcd*gcd(apart_co, power(p_gcd, exp_b-1)*power(p_co, exp_b), ca, cb, false);
+ } // p_gcd.is_equal(_ex1)
}
#endif