Source code for sympy.polys.euclidtools

"""Euclidean algorithms, GCDs, LCMs and polynomial remainder sequences. """

from sympy.polys.densebasic import (
    dup_strip, dmp_raise,
    dmp_zero, dmp_one, dmp_ground,
    dmp_one_p, dmp_zero_p,
    dmp_zeros,
    dup_degree, dmp_degree, dmp_degree_in,
    dup_LC, dmp_LC, dmp_ground_LC,
    dmp_multi_deflate, dmp_inflate,
    dup_convert, dmp_convert,
    dmp_apply_pairs)

from sympy.polys.densearith import (
    dup_sub_mul,
    dup_neg, dmp_neg,
    dmp_add,
    dmp_sub,
    dup_mul, dmp_mul,
    dmp_pow,
    dup_div, dmp_div,
    dup_rem,
    dup_quo, dmp_quo,
    dup_prem, dmp_prem,
    dup_mul_ground, dmp_mul_ground,
    dmp_mul_term,
    dup_quo_ground, dmp_quo_ground,
    dup_max_norm, dmp_max_norm)

from sympy.polys.densetools import (
    dup_clear_denoms, dmp_clear_denoms,
    dup_diff, dmp_diff,
    dup_eval, dmp_eval, dmp_eval_in,
    dup_trunc, dmp_ground_trunc,
    dup_monic, dmp_ground_monic,
    dup_primitive, dmp_ground_primitive,
    dup_extract, dmp_ground_extract)

from sympy.polys.galoistools import (
    gf_int, gf_crt)

from sympy.polys.polyerrors import (
    MultivariatePolynomialError,
    HeuristicGCDFailed,
    HomomorphismFailed,
    NotInvertible,
    DomainError)

from sympy.polys.polyconfig import query

from sympy.utilities import cythonized

from sympy.ntheory import nextprime

def dup_half_gcdex(f, g, K):
    """
    Half extended Euclidean algorithm in `F[x]`.

    Returns ``(s, h)`` such that ``h = gcd(f, g)`` and ``s*f = h (mod g)``.

    Examples
    ========

    >>> from sympy.polys.domains import QQ
    >>> from sympy.polys.euclidtools import dup_half_gcdex

    >>> f = QQ.map([1, -2, -6, 12, 15])
    >>> g = QQ.map([1, 1, -4, -4])

    >>> dup_half_gcdex(f, g, QQ)
    ([-1/5, 3/5], [1/1, 1/1])

    """
    if not (K.has_Field or not K.is_Exact):
        raise DomainError("can't compute half extended GCD over %s" % K)

    a, b = [K.one], []

    while g:
        q, r = dup_div(f, g, K)
        f, g = g, r
        a, b = b, dup_sub_mul(a, q, b, K)

    a = dup_quo_ground(a, dup_LC(f, K), K)
    f = dup_monic(f, K)

    return a, f

[docs]def dmp_half_gcdex(f, g, u, K): """ Half extended Euclidean algorithm in `F[X]`. Examples ======== >>> from sympy.polys.domains import QQ >>> from sympy.polys.euclidtools import dmp_half_gcdex """ if not u: return dup_half_gcdex(f, g, K) else: raise MultivariatePolynomialError(f, g)
def dup_gcdex(f, g, K): """ Extended Euclidean algorithm in `F[x]`. Returns ``(s, t, h)`` such that ``h = gcd(f, g)`` and ``s*f + t*g = h``. Examples ======== >>> from sympy.polys.domains import QQ >>> from sympy.polys.euclidtools import dup_gcdex >>> f = QQ.map([1, -2, -6, 12, 15]) >>> g = QQ.map([1, 1, -4, -4]) >>> dup_gcdex(f, g, QQ) ([-1/5, 3/5], [1/5, -6/5, 2/1], [1/1, 1/1]) """ s, h = dup_half_gcdex(f, g, K) F = dup_sub_mul(h, s, f, K) t = dup_quo(F, g, K) return s, t, h
[docs]def dmp_gcdex(f, g, u, K): """ Extended Euclidean algorithm in `F[X]`. Examples ======== >>> from sympy.polys.domains import QQ >>> from sympy.polys.euclidtools import dmp_gcdex """ if not u: return dup_gcdex(f, g, K) else: raise MultivariatePolynomialError(f, g)
def dup_invert(f, g, K): """ Compute multiplicative inverse of `f` modulo `g` in `F[x]`. Examples ======== >>> from sympy.polys.domains import QQ >>> from sympy.polys.euclidtools import dup_invert >>> f = QQ.map([1, 0, -1]) >>> g = QQ.map([2, -1]) >>> h = QQ.map([1, -1]) >>> dup_invert(f, g, QQ) [-4/3] >>> dup_invert(f, h, QQ) Traceback (most recent call last): ... NotInvertible: zero divisor """ s, h = dup_half_gcdex(f, g, K) if h == [K.one]: return dup_rem(s, g, K) else: raise NotInvertible("zero divisor")
[docs]def dmp_invert(f, g, u, K): """ Compute multiplicative inverse of `f` modulo `g` in `F[X]`. Examples ======== >>> from sympy.polys.domains import QQ >>> from sympy.polys.euclidtools import dmp_invert """ if not u: return dup_invert(f, g, K) else: raise MultivariatePolynomialError(f, g)
def dup_euclidean_prs(f, g, K): """ Euclidean polynomial remainder sequence (PRS) in `K[x]`. Examples ======== >>> from sympy.polys.domains import QQ >>> from sympy.polys.euclidtools import dup_euclidean_prs >>> f = QQ.map([1, 0, 1, 0, -3, -3, 8, 2, -5]) >>> g = QQ.map([3, 0, 5, 0, -4, -9, 21]) >>> prs = dup_euclidean_prs(f, g, QQ) >>> prs[0] [1/1, 0/1, 1/1, 0/1, -3/1, -3/1, 8/1, 2/1, -5/1] >>> prs[1] [3/1, 0/1, 5/1, 0/1, -4/1, -9/1, 21/1] >>> prs[2] [-5/9, 0/1, 1/9, 0/1, -1/3] >>> prs[3] [-117/25, -9/1, 441/25] >>> prs[4] [233150/19773, -102500/6591] >>> prs[5] [-1288744821/543589225] """ prs = [f, g] h = dup_rem(f, g, K) while h: prs.append(h) f, g = g, h h = dup_rem(f, g, K) return prs
[docs]def dmp_euclidean_prs(f, g, u, K): """ Euclidean polynomial remainder sequence (PRS) in `K[X]`. Examples ======== >>> from sympy.polys.domains import QQ >>> from sympy.polys.euclidtools import dmp_euclidean_prs """ if not u: return dup_euclidean_prs(f, g, K) else: raise MultivariatePolynomialError(f, g)
def dup_primitive_prs(f, g, K): """ Primitive polynomial remainder sequence (PRS) in `K[x]`. Examples ======== >>> from sympy.polys.domains import ZZ >>> from sympy.polys.euclidtools import dup_primitive_prs >>> f = ZZ.map([1, 0, 1, 0, -3, -3, 8, 2, -5]) >>> g = ZZ.map([3, 0, 5, 0, -4, -9, 21]) >>> prs = dup_primitive_prs(f, g, ZZ) >>> prs[0] [1, 0, 1, 0, -3, -3, 8, 2, -5] >>> prs[1] [3, 0, 5, 0, -4, -9, 21] >>> prs[2] [-5, 0, 1, 0, -3] >>> prs[3] [13, 25, -49] >>> prs[4] [4663, -6150] >>> prs[5] [1] """ prs = [f, g] _, h = dup_primitive(dup_prem(f, g, K), K) while h: prs.append(h) f, g = g, h _, h = dup_primitive(dup_prem(f, g, K), K) return prs
[docs]def dmp_primitive_prs(f, g, u, K): """ Primitive polynomial remainder sequence (PRS) in `K[X]`. Examples ======== >>> from sympy.polys.domains import ZZ >>> from sympy.polys.euclidtools import dmp_primitive_prs """ if not u: return dup_primitive_prs(f, g, K) else: raise MultivariatePolynomialError(f, g)
@cythonized("n,m,d,k") def dup_inner_subresultants(f, g, K): """ Subresultant PRS algorithm in `K[x]`. Computes the subresultant polynomial remainder sequence (PRS) of `f` and `g`, and the values for `\beta_i` and `\delta_i`. The last two sequences of values are necessary for computing the resultant in :func:`dup_prs_resultant`. Examples ======== >>> from sympy.polys.domains import ZZ >>> from sympy.polys.euclidtools import dup_inner_subresultants >>> f = ZZ.map([1, 0, 1]) >>> g = ZZ.map([1, 0, -1]) >>> dup_inner_subresultants(f, g, ZZ) ([[1, 0, 1], [1, 0, -1], [-2]], [-1, -1], [0, 2]) """ n = dup_degree(f) m = dup_degree(g) if n < m: f, g = g, f n, m = m, n R = [f, g] d = n - m b = (-K.one)**(d+1) c = -K.one B, D = [b], [d] if not f or not g: return R, B, D h = dup_prem(f, g, K) h = dup_mul_ground(h, b, K) while h: k = dup_degree(h) R.append(h) lc = dup_LC(g, K) if not d: q = c else: q = c**(d-1) c = K.quo((-lc)**d, q) b = -lc * c**(m-k) f, g, m, d = g, h, k, m-k B.append(b) D.append(d) h = dup_prem(f, g, K) h = dup_quo_ground(h, b, K) return R, B, D def dup_subresultants(f, g, K): """ Computes subresultant PRS of two polynomials in `K[x]`. Examples ======== >>> from sympy.polys.domains import ZZ >>> from sympy.polys.euclidtools import dup_subresultants >>> f = ZZ.map([1, 0, 1]) >>> g = ZZ.map([1, 0, -1]) >>> dup_subresultants(f, g, ZZ) [[1, 0, 1], [1, 0, -1], [-2]] """ return dup_inner_subresultants(f, g, K)[0] @cythonized("s,i,du,dv,dw") def dup_prs_resultant(f, g, K): """ Resultant algorithm in `K[x]` using subresultant PRS. Examples ======== >>> from sympy.polys.domains import ZZ >>> from sympy.polys.euclidtools import dup_prs_resultant >>> f = ZZ.map([1, 0, 1]) >>> g = ZZ.map([1, 0, -1]) >>> dup_prs_resultant(f, g, ZZ) (4, [[1, 0, 1], [1, 0, -1], [-2]]) """ if not f or not g: return (K.zero, []) R, B, D = dup_inner_subresultants(f, g, K) if dup_degree(R[-1]) > 0: return (K.zero, R) if R[-2] == [K.one]: return (dup_LC(R[-1], K), R) s, i = 1, 1 p, q = K.one, K.one for b, d in list(zip(B, D))[:-1]: du = dup_degree(R[i-1]) dv = dup_degree(R[i ]) dw = dup_degree(R[i+1]) if du % 2 and dv % 2: s = -s lc, i = dup_LC(R[i], K), i+1 p *= b**dv * lc**(du-dw) q *= lc**(dv*(1+d)) if s < 0: p = -p i = dup_degree(R[-2]) res = dup_LC(R[-1], K)**i res = K.quo(res*p, q) return res, R def dup_resultant(f, g, K): """ Computes resultant of two polynomials in `K[x]`. Examples ======== >>> from sympy.polys.domains import ZZ >>> from sympy.polys.euclidtools import dup_resultant >>> f = ZZ.map([1, 0, 1]) >>> g = ZZ.map([1, 0, -1]) >>> dup_resultant(f, g, ZZ) 4 """ return dup_prs_resultant(f, g, K)[0] @cythonized("u,v,n,m,d,k")
[docs]def dmp_inner_subresultants(f, g, u, K): """ Subresultant PRS algorithm in `K[X]`. Examples ======== >>> from sympy.polys.domains import ZZ >>> from sympy.polys.euclidtools import dmp_inner_subresultants >>> f = ZZ.map([[3, 0], [], [-1, 0, 0, -4]]) >>> g = ZZ.map([[1], [1, 0, 0, 0], [-9]]) >>> a = [[3, 0, 0, 0, 0], [1, 0, -27, 4]] >>> b = [[-3, 0, 0, -12, 1, 0, -54, 8, 729, -216, 16]] >>> R = ZZ.map([f, g, a, b]) >>> B = ZZ.map([[-1], [1], [9, 0, 0, 0, 0, 0, 0, 0, 0]]) >>> D = ZZ.map([0, 1, 1]) >>> dmp_inner_subresultants(f, g, 1, ZZ) == (R, B, D) True """ if not u: return dup_inner_subresultants(f, g, K) n = dmp_degree(f, u) m = dmp_degree(g, u) if n < m: f, g = g, f n, m = m, n R = [f, g] d = n - m v = u - 1 b = dmp_pow(dmp_ground(-K.one, v), d+1, v, K) c = dmp_ground(-K.one, v) B, D = [b], [d] if dmp_zero_p(f, u) or dmp_zero_p(g, u): return R, B, D h = dmp_prem(f, g, u, K) h = dmp_mul_term(h, b, 0, u, K) while not dmp_zero_p(h, u): k = dmp_degree(h, u) R.append(h) lc = dmp_LC(g, K) p = dmp_pow(dmp_neg(lc, v, K), d, v, K) if not d: q = c else: q = dmp_pow(c, d-1, v, K) c = dmp_quo(p, q, v, K) b = dmp_mul(dmp_neg(lc, v, K), dmp_pow(c, m-k, v, K), v, K) f, g, m, d = g, h, k, m-k B.append(b) D.append(d) h = dmp_prem(f, g, u, K) h = [ dmp_quo(ch, b, v, K) for ch in h ] return R, B, D
@cythonized("u")
[docs]def dmp_subresultants(f, g, u, K): """ Computes subresultant PRS of two polynomials in `K[X]`. Examples ======== >>> from sympy.polys.domains import ZZ >>> from sympy.polys.euclidtools import dmp_subresultants >>> f = [[3, 0], [], [-1, 0, 0, -4]] >>> g = [[1], [1, 0, 0, 0], [-9]] >>> a = [[3, 0, 0, 0, 0], [1, 0, -27, 4]] >>> b = [[-3, 0, 0, -12, 1, 0, -54, 8, 729, -216, 16]] >>> dmp_subresultants(f, g, 1, ZZ) == [f, g, a, b] True """ return dmp_inner_subresultants(f, g, u, K)[0]
@cythonized("u,v,s,i,d,du,dv,dw")
[docs]def dmp_prs_resultant(f, g, u, K): """ Resultant algorithm in `K[X]` using subresultant PRS. Examples ======== >>> from sympy.polys.domains import ZZ >>> from sympy.polys.euclidtools import dmp_prs_resultant >>> f = ZZ.map([[3, 0], [], [-1, 0, 0, -4]]) >>> g = ZZ.map([[1], [1, 0, 0, 0], [-9]]) >>> a = ZZ.map([[3, 0, 0, 0, 0], [1, 0, -27, 4]]) >>> b = ZZ.map([[-3, 0, 0, -12, 1, 0, -54, 8, 729, -216, 16]]) >>> dmp_prs_resultant(f, g, 1, ZZ) == (b[0], [f, g, a, b]) True """ if not u: return dup_prs_resultant(f, g, K) if dmp_zero_p(f, u) or dmp_zero_p(g, u): return (dmp_zero(u-1), []) R, B, D = dmp_inner_subresultants(f, g, u, K) if dmp_degree(R[-1], u) > 0: return (dmp_zero(u-1), R) if dmp_one_p(R[-2], u, K): return (dmp_LC(R[-1], K), R) s, i, v = 1, 1, u-1 p = dmp_one(v, K) q = dmp_one(v, K) for b, d in list(zip(B, D))[:-1]: du = dmp_degree(R[i-1], u) dv = dmp_degree(R[i ], u) dw = dmp_degree(R[i+1], u) if du % 2 and dv % 2: s = -s lc, i = dmp_LC(R[i], K), i+1 p = dmp_mul(dmp_mul(p, dmp_pow(b, dv, v, K), v, K), dmp_pow(lc, du-dw, v, K), v, K) q = dmp_mul(q, dmp_pow(lc, dv*(1+d), v, K), v, K) _, p, q = dmp_inner_gcd(p, q, v, K) if s < 0: p = dmp_neg(p, v, K) i = dmp_degree(R[-2], u) res = dmp_pow(dmp_LC(R[-1], K), i, v, K) res = dmp_quo(dmp_mul(res, p, v, K), q, v, K) return res, R
@cythonized("u,v,n,m,N,M,B")
[docs]def dmp_zz_modular_resultant(f, g, p, u, K): """ Compute resultant of `f` and `g` modulo a prime `p`. Examples ======== >>> from sympy.polys.domains import ZZ >>> from sympy.polys.euclidtools import dmp_zz_modular_resultant >>> f = ZZ.map([[1], [1, 2]]) >>> g = ZZ.map([[2, 1], [3]]) >>> dmp_zz_modular_resultant(f, g, ZZ(5), 1, ZZ) [-2, 0, 1] """ if not u: return gf_int(dup_prs_resultant(f, g, K)[0] % p, p) v = u - 1 n = dmp_degree(f, u) m = dmp_degree(g, u) N = dmp_degree_in(f, 1, u) M = dmp_degree_in(g, 1, u) B = n*M + m*N D, a = [K.one], -K.one r = dmp_zero(v) while dup_degree(D) <= B: while True: a += K.one if a == p: raise HomomorphismFailed('no luck') F = dmp_eval_in(f, gf_int(a, p), 1, u, K) if dmp_degree(F, v) == n: G = dmp_eval_in(g, gf_int(a, p), 1, u, K) if dmp_degree(G, v) == m: break R = dmp_zz_modular_resultant(F, G, p, v, K) e = dmp_eval(r, a, v, K) if not v: R = dup_strip([R]) e = dup_strip([e]) else: R = [R] e = [e] d = K.invert(dup_eval(D, a, K), p) d = dup_mul_ground(D, d, K) d = dmp_raise(d, v, 0, K) c = dmp_mul(d, dmp_sub(R, e, v, K), v, K) r = dmp_add(r, c, v, K) r = dmp_ground_trunc(r, p, v, K) D = dup_mul(D, [K.one, -a], K) D = dup_trunc(D, p, K) return r
def _collins_crt(r, R, P, p, K): """Wrapper of CRT for Collins's resultant algorithm. """ return gf_int(gf_crt([r, R], [P, p], K), P*p) @cythonized("u,v,n,m")
[docs]def dmp_zz_collins_resultant(f, g, u, K): """ Collins's modular resultant algorithm in `Z[X]`. Examples ======== >>> from sympy.polys.domains import ZZ >>> from sympy.polys.euclidtools import dmp_zz_collins_resultant >>> f = ZZ.map([[1], [1, 2]]) >>> g = ZZ.map([[2, 1], [3]]) >>> dmp_zz_collins_resultant(f, g, 1, ZZ) [-2, -5, 1] """ n = dmp_degree(f, u) m = dmp_degree(g, u) if n < 0 or m < 0: return dmp_zero(u-1) A = dmp_max_norm(f, u, K) B = dmp_max_norm(g, u, K) a = dmp_ground_LC(f, u, K) b = dmp_ground_LC(g, u, K) v = u - 1 B = K(2)*K.factorial(n+m)*A**m*B**n r, p, P = dmp_zero(v), K.one, K.one while P <= B: p = K(nextprime(p)) while not (a % p) or not (b % p): p = K(nextprime(p)) F = dmp_ground_trunc(f, p, u, K) G = dmp_ground_trunc(g, p, u, K) try: R = dmp_zz_modular_resultant(F, G, p, u, K) except HomomorphismFailed: continue if K.is_one(P): r = R else: r = dmp_apply_pairs(r, R, _collins_crt, (P, p, K), v, K) P *= p return r
@cythonized("u,n,m")
[docs]def dmp_qq_collins_resultant(f, g, u, K0): """ Collins's modular resultant algorithm in `Q[X]`. Examples ======== >>> from sympy.polys.domains import QQ >>> from sympy.polys.euclidtools import dmp_qq_collins_resultant >>> f = [[QQ(1,2)], [QQ(1), QQ(2,3)]] >>> g = [[QQ(2), QQ(1)], [QQ(3)]] >>> dmp_qq_collins_resultant(f, g, 1, QQ) [-2/1, -7/3, 5/6] """ n = dmp_degree(f, u) m = dmp_degree(g, u) if n < 0 or m < 0: return dmp_zero(u-1) K1 = K0.get_ring() cf, f = dmp_clear_denoms(f, u, K0, K1) cg, g = dmp_clear_denoms(g, u, K0, K1) f = dmp_convert(f, u, K0, K1) g = dmp_convert(g, u, K0, K1) r = dmp_zz_collins_resultant(f, g, u, K1) r = dmp_convert(r, u-1, K1, K0) c = K0.convert(cf**m * cg**n, K1) return dmp_quo_ground(r, c, u-1, K0)
@cythonized("u")
[docs]def dmp_resultant(f, g, u, K): """ Computes resultant of two polynomials in `K[X]`. Examples ======== >>> from sympy.polys.domains import ZZ >>> from sympy.polys.euclidtools import dmp_resultant >>> f = ZZ.map([[3, 0], [], [-1, 0, 0, -4]]) >>> g = ZZ.map([[1], [1, 0, 0, 0], [-9]]) >>> dmp_resultant(f, g, 1, ZZ) [-3, 0, 0, -12, 1, 0, -54, 8, 729, -216, 16] """ if not u: return dup_resultant(f, g, K) if K.has_Field: if K.is_QQ and query('USE_COLLINS_RESULTANT'): return dmp_qq_collins_resultant(f, g, u, K) else: if K.is_ZZ and query('USE_COLLINS_RESULTANT'): return dmp_zz_collins_resultant(f, g, u, K) return dmp_prs_resultant(f, g, u, K)[0]
@cythonized("d,s") def dup_discriminant(f, K): """ Computes discriminant of a polynomial in `K[x]`. Examples ======== >>> from sympy.polys.domains import ZZ >>> from sympy.polys.euclidtools import dup_discriminant >>> dup_discriminant([ZZ(1), ZZ(2), ZZ(3)], ZZ) -8 """ d = dup_degree(f) if d <= 0: return K.zero else: s = (-1)**((d*(d-1)) // 2) c = dup_LC(f, K) r = dup_resultant(f, dup_diff(f, 1, K), K) return K.quo(r, c*K(s)) @cythonized("u,v,d,s")
[docs]def dmp_discriminant(f, u, K): """ Computes discriminant of a polynomial in `K[X]`. Examples ======== >>> from sympy.polys.domains import ZZ >>> from sympy.polys.euclidtools import dmp_discriminant >>> f = ZZ.map([[[[1]], [[]]], [[[1], []]], [[[1, 0]]]]) >>> dmp_discriminant(f, 3, ZZ) [[[-4, 0]], [[1], [], []]] """ if not u: return dup_discriminant(f, K) d, v = dmp_degree(f, u), u-1 if d <= 0: return dmp_zero(v) else: s = (-1)**((d*(d-1)) // 2) c = dmp_LC(f, K) r = dmp_resultant(f, dmp_diff(f, 1, u, K), u, K) c = dmp_mul_ground(c, K(s), v, K) return dmp_quo(r, c, v, K)
def _dup_rr_trivial_gcd(f, g, K): """Handle trivial cases in GCD algorithm over a ring. """ if not (f or g): return [], [], [] elif not f: if K.is_nonnegative(dup_LC(g, K)): return g, [], [K.one] else: return dup_neg(g, K), [], [-K.one] elif not g: if K.is_nonnegative(dup_LC(f, K)): return f, [K.one], [] else: return dup_neg(f, K), [-K.one], [] return None def _dup_ff_trivial_gcd(f, g, K): """Handle trivial cases in GCD algorithm over a field. """ if not (f or g): return [], [], [] elif not f: return dup_monic(g, K), [], [dup_LC(g, K)] elif not g: return dup_monic(f, K), [dup_LC(f, K)], [] else: return None @cythonized("u") def _dmp_rr_trivial_gcd(f, g, u, K): """Handle trivial cases in GCD algorithm over a ring. """ zero_f = dmp_zero_p(f, u) zero_g = dmp_zero_p(g, u) if zero_f and zero_g: return tuple(dmp_zeros(3, u, K)) elif zero_f: if K.is_nonnegative(dmp_ground_LC(g, u, K)): return g, dmp_zero(u), dmp_one(u, K) else: return dmp_neg(g, u, K), dmp_zero(u), dmp_ground(-K.one, u) elif zero_g: if K.is_nonnegative(dmp_ground_LC(f, u, K)): return f, dmp_one(u, K), dmp_zero(u) else: return dmp_neg(f, u, K), dmp_ground(-K.one, u), dmp_zero(u) elif query('USE_SIMPLIFY_GCD'): return _dmp_simplify_gcd(f, g, u, K) else: return None @cythonized("u") def _dmp_ff_trivial_gcd(f, g, u, K): """Handle trivial cases in GCD algorithm over a field. """ zero_f = dmp_zero_p(f, u) zero_g = dmp_zero_p(g, u) if zero_f and zero_g: return tuple(dmp_zeros(3, u, K)) elif zero_f: return (dmp_ground_monic(g, u, K), dmp_zero(u), dmp_ground(dmp_ground_LC(g, u, K), u)) elif zero_g: return (dmp_ground_monic(f, u, K), dmp_ground(dmp_ground_LC(f, u, K), u), dmp_zero(u)) elif query('USE_SIMPLIFY_GCD'): return _dmp_simplify_gcd(f, g, u, K) else: return None @cythonized("u,v,df,dg") def _dmp_simplify_gcd(f, g, u, K): """Try to eliminate `x_0` from GCD computation in `K[X]`. """ df = dmp_degree(f, u) dg = dmp_degree(g, u) if df > 0 and dg > 0: return None if not (df or dg): F = dmp_LC(f, K) G = dmp_LC(g, K) else: if not df: F = dmp_LC(f, K) G = dmp_content(g, u, K) else: F = dmp_content(f, u, K) G = dmp_LC(g, K) v = u - 1 h = dmp_gcd(F, G, v, K) cff = [ dmp_quo(cf, h, v, K) for cf in f ] cfg = [ dmp_quo(cg, h, v, K) for cg in g ] return [h], cff, cfg def dup_rr_prs_gcd(f, g, K): """ Computes polynomial GCD using subresultants over a ring. Returns ``(h, cff, cfg)`` such that ``a = gcd(f, g)``, ``cff = quo(f, h)``, and ``cfg = quo(g, h)``. Examples ======== >>> from sympy.polys.domains import ZZ >>> from sympy.polys.euclidtools import dup_rr_prs_gcd >>> f = ZZ.map([1, 0, -1]) >>> g = ZZ.map([1, -3, 2]) >>> dup_rr_prs_gcd(f, g, ZZ) ([1, -1], [1, 1], [1, -2]) """ result = _dup_rr_trivial_gcd(f, g, K) if result is not None: return result fc, F = dup_primitive(f, K) gc, G = dup_primitive(g, K) c = K.gcd(fc, gc) h = dup_subresultants(F, G, K)[-1] _, h = dup_primitive(h, K) if K.is_negative(dup_LC(h, K)): c = -c h = dup_mul_ground(h, c, K) cff = dup_quo(f, h, K) cfg = dup_quo(g, h, K) return h, cff, cfg def dup_ff_prs_gcd(f, g, K): """ Computes polynomial GCD using subresultants over a field. Returns ``(h, cff, cfg)`` such that ``a = gcd(f, g)``, ``cff = quo(f, h)``, and ``cfg = quo(g, h)``. Examples ======== >>> from sympy.polys.domains import QQ >>> from sympy.polys.euclidtools import dup_ff_prs_gcd >>> f = QQ.map([1, 0, -1]) >>> g = QQ.map([1, -3, 2]) >>> dup_ff_prs_gcd(f, g, QQ) ([1/1, -1/1], [1/1, 1/1], [1/1, -2/1]) """ result = _dup_ff_trivial_gcd(f, g, K) if result is not None: return result h = dup_subresultants(f, g, K)[-1] h = dup_monic(h, K) cff = dup_quo(f, h, K) cfg = dup_quo(g, h, K) return h, cff, cfg @cythonized("u")
[docs]def dmp_rr_prs_gcd(f, g, u, K): """ Computes polynomial GCD using subresultants over a ring. Returns ``(h, cff, cfg)`` such that ``a = gcd(f, g)``, ``cff = quo(f, h)``, and ``cfg = quo(g, h)``. Examples ======== >>> from sympy.polys.domains import ZZ >>> from sympy.polys.euclidtools import dmp_rr_prs_gcd >>> f = ZZ.map([[1], [2, 0], [1, 0, 0]]) >>> g = ZZ.map([[1], [1, 0], []]) >>> dmp_rr_prs_gcd(f, g, 1, ZZ) ([[1], [1, 0]], [[1], [1, 0]], [[1], []]) """ if not u: return dup_rr_prs_gcd(f, g, K) result = _dmp_rr_trivial_gcd(f, g, u, K) if result is not None: return result fc, F = dmp_primitive(f, u, K) gc, G = dmp_primitive(g, u, K) h = dmp_subresultants(F, G, u, K)[-1] c, _, _ = dmp_rr_prs_gcd(fc, gc, u-1, K) if K.is_negative(dmp_ground_LC(h, u, K)): h = dmp_neg(h, u, K) _, h = dmp_primitive(h, u, K) h = dmp_mul_term(h, c, 0, u, K) cff = dmp_quo(f, h, u, K) cfg = dmp_quo(g, h, u, K) return h, cff, cfg
@cythonized("u")
[docs]def dmp_ff_prs_gcd(f, g, u, K): """ Computes polynomial GCD using subresultants over a field. Returns ``(h, cff, cfg)`` such that ``a = gcd(f, g)``, ``cff = quo(f, h)``, and ``cfg = quo(g, h)``. Examples ======== >>> from sympy.polys.domains import QQ >>> from sympy.polys.euclidtools import dmp_ff_prs_gcd >>> f = [[QQ(1,2)], [QQ(1), QQ(0)], [QQ(1,2), QQ(0), QQ(0)]] >>> g = [[QQ(1)], [QQ(1), QQ(0)], []] >>> dmp_ff_prs_gcd(f, g, 1, QQ) ([[1/1], [1/1, 0/1]], [[1/2], [1/2, 0/1]], [[1/1], []]) """ if not u: return dup_ff_prs_gcd(f, g, K) result = _dmp_ff_trivial_gcd(f, g, u, K) if result is not None: return result fc, F = dmp_primitive(f, u, K) gc, G = dmp_primitive(g, u, K) h = dmp_subresultants(F, G, u, K)[-1] c, _, _ = dmp_ff_prs_gcd(fc, gc, u-1, K) _, h = dmp_primitive(h, u, K) h = dmp_mul_term(h, c, 0, u, K) h = dmp_ground_monic(h, u, K) cff = dmp_quo(f, h, u, K) cfg = dmp_quo(g, h, u, K) return h, cff, cfg
HEU_GCD_MAX = 6 def _dup_zz_gcd_interpolate(h, x, K): """Interpolate polynomial GCD from integer GCD. """ f = [] while h: g = h % x if g > x // 2: g -= x f.insert(0, g) h = (h-g) // x return f @cythonized("i,df,dg") def dup_zz_heu_gcd(f, g, K): """ Heuristic polynomial GCD in `Z[x]`. Given univariate polynomials `f` and `g` in `Z[x]`, returns their GCD and cofactors, i.e. polynomials ``h``, ``cff`` and ``cfg`` such that:: h = gcd(f, g), cff = quo(f, h) and cfg = quo(g, h) The algorithm is purely heuristic which means it may fail to compute the GCD. This will be signaled by raising an exception. In this case you will need to switch to another GCD method. The algorithm computes the polynomial GCD by evaluating polynomials f and g at certain points and computing (fast) integer GCD of those evaluations. The polynomial GCD is recovered from the integer image by interpolation. The final step is to verify if the result is the correct GCD. This gives cofactors as a side effect. Examples ======== >>> from sympy.polys.domains import ZZ >>> from sympy.polys.euclidtools import dup_zz_heu_gcd >>> f = ZZ.map([1, 0, -1]) >>> g = ZZ.map([1, -3, 2]) >>> dup_zz_heu_gcd(f, g, ZZ) ([1, -1], [1, 1], [1, -2]) References ========== 1. [Liao95]_ """ result = _dup_rr_trivial_gcd(f, g, K) if result is not None: return result df = dup_degree(f) dg = dup_degree(g) gcd, f, g = dup_extract(f, g, K) if df == 0 or dg == 0: return [gcd], f, g f_norm = dup_max_norm(f, K) g_norm = dup_max_norm(g, K) B = 2*min(f_norm, g_norm) + 29 x = max(min(B, 99*K.sqrt(B)), 2*min(f_norm // abs(dup_LC(f, K)), g_norm // abs(dup_LC(g, K))) + 2) for i in xrange(0, HEU_GCD_MAX): ff = dup_eval(f, x, K) gg = dup_eval(g, x, K) if ff and gg: h = K.gcd(ff, gg) cff = ff // h cfg = gg // h h = _dup_zz_gcd_interpolate(h, x, K) h = dup_primitive(h, K)[1] cff_, r = dup_div(f, h, K) if not r: cfg_, r = dup_div(g, h, K) if not r: h = dup_mul_ground(h, gcd, K) return h, cff_, cfg_ cff = _dup_zz_gcd_interpolate(cff, x, K) h, r = dup_div(f, cff, K) if not r: cfg_, r = dup_div(g, h, K) if not r: h = dup_mul_ground(h, gcd, K) return h, cff, cfg_ cfg = _dup_zz_gcd_interpolate(cfg, x, K) h, r = dup_div(g, cfg, K) if not r: cff_, r = dup_div(f, h, K) if not r: h = dup_mul_ground(h, gcd, K) return h, cff_, cfg x = 73794*x * K.sqrt(K.sqrt(x)) // 27011 raise HeuristicGCDFailed('no luck') @cythonized("v") def _dmp_zz_gcd_interpolate(h, x, v, K): """Interpolate polynomial GCD from integer GCD. """ f = [] while not dmp_zero_p(h, v): g = dmp_ground_trunc(h, x, v, K) f.insert(0, g) h = dmp_sub(h, g, v, K) h = dmp_quo_ground(h, x, v, K) if K.is_negative(dmp_ground_LC(f, v+1, K)): return dmp_neg(f, v+1, K) else: return f @cythonized("u,v,i,dg,df")
[docs]def dmp_zz_heu_gcd(f, g, u, K): """ Heuristic polynomial GCD in `Z[X]`. Given univariate polynomials `f` and `g` in `Z[X]`, returns their GCD and cofactors, i.e. polynomials ``h``, ``cff`` and ``cfg`` such that:: h = gcd(f, g), cff = quo(f, h) and cfg = quo(g, h) The algorithm is purely heuristic which means it may fail to compute the GCD. This will be signaled by raising an exception. In this case you will need to switch to another GCD method. The algorithm computes the polynomial GCD by evaluating polynomials f and g at certain points and computing (fast) integer GCD of those evaluations. The polynomial GCD is recovered from the integer image by interpolation. The evaluation proces reduces f and g variable by variable into a large integer. The final step is to verify if the interpolated polynomial is the correct GCD. This gives cofactors of the input polynomials as a side effect. Examples ======== >>> from sympy.polys.domains import ZZ >>> from sympy.polys.euclidtools import dmp_zz_heu_gcd >>> f = ZZ.map([[1], [2, 0], [1, 0, 0]]) >>> g = ZZ.map([[1], [1, 0], []]) >>> dmp_zz_heu_gcd(f, g, 1, ZZ) ([[1], [1, 0]], [[1], [1, 0]], [[1], []]) References ========== 1. [Liao95]_ """ if not u: return dup_zz_heu_gcd(f, g, K) result = _dmp_rr_trivial_gcd(f, g, u, K) if result is not None: return result gcd, f, g = dmp_ground_extract(f, g, u, K) f_norm = dmp_max_norm(f, u, K) g_norm = dmp_max_norm(g, u, K) B = 2*min(f_norm, g_norm) + 29 x = max(min(B, 99*K.sqrt(B)), 2*min(f_norm // abs(dmp_ground_LC(f, u, K)), g_norm // abs(dmp_ground_LC(g, u, K))) + 2) for i in xrange(0, HEU_GCD_MAX): ff = dmp_eval(f, x, u, K) gg = dmp_eval(g, x, u, K) v = u - 1 if not (dmp_zero_p(ff, v) or dmp_zero_p(gg, v)): h, cff, cfg = dmp_zz_heu_gcd(ff, gg, v, K) h = _dmp_zz_gcd_interpolate(h, x, v, K) h = dmp_ground_primitive(h, u, K)[1] cff_, r = dmp_div(f, h, u, K) if dmp_zero_p(r, u): cfg_, r = dmp_div(g, h, u, K) if dmp_zero_p(r, u): h = dmp_mul_ground(h, gcd, u, K) return h, cff_, cfg_ cff = _dmp_zz_gcd_interpolate(cff, x, v, K) h, r = dmp_div(f, cff, u, K) if dmp_zero_p(r, u): cfg_, r = dmp_div(g, h, u, K) if dmp_zero_p(r, u): h = dmp_mul_ground(h, gcd, u, K) return h, cff, cfg_ cfg = _dmp_zz_gcd_interpolate(cfg, x, v, K) h, r = dmp_div(g, cfg, u, K) if dmp_zero_p(r, u): cff_, r = dmp_div(f, h, u, K) if dmp_zero_p(r, u): h = dmp_mul_ground(h, gcd, u, K) return h, cff_, cfg x = 73794*x * K.sqrt(K.sqrt(x)) // 27011 raise HeuristicGCDFailed('no luck')
def dup_qq_heu_gcd(f, g, K0): """ Heuristic polynomial GCD in `Q[x]`. Returns ``(h, cff, cfg)`` such that ``a = gcd(f, g)``, ``cff = quo(f, h)``, and ``cfg = quo(g, h)``. Examples ======== >>> from sympy.polys.domains import QQ >>> from sympy.polys.euclidtools import dup_qq_heu_gcd >>> f = [QQ(1,2), QQ(7,4), QQ(3,2)] >>> g = [QQ(1,2), QQ(1), QQ(0)] >>> dup_qq_heu_gcd(f, g, QQ) ([1/1, 2/1], [1/2, 3/4], [1/2, 0/1]) """ result = _dup_ff_trivial_gcd(f, g, K0) if result is not None: return result K1 = K0.get_ring() cf, f = dup_clear_denoms(f, K0, K1) cg, g = dup_clear_denoms(g, K0, K1) f = dup_convert(f, K0, K1) g = dup_convert(g, K0, K1) h, cff, cfg = dup_zz_heu_gcd(f, g, K1) h = dup_convert(h, K1, K0) c = dup_LC(h, K0) h = dup_monic(h, K0) cff = dup_convert(cff, K1, K0) cfg = dup_convert(cfg, K1, K0) cff = dup_mul_ground(cff, K0.quo(c, cf), K0) cfg = dup_mul_ground(cfg, K0.quo(c, cg), K0) return h, cff, cfg @cythonized("u")
[docs]def dmp_qq_heu_gcd(f, g, u, K0): """ Heuristic polynomial GCD in `Q[X]`. Returns ``(h, cff, cfg)`` such that ``a = gcd(f, g)``, ``cff = quo(f, h)``, and ``cfg = quo(g, h)``. Examples ======== >>> from sympy.polys.domains import QQ >>> from sympy.polys.euclidtools import dmp_qq_heu_gcd >>> f = [[QQ(1,4)], [QQ(1), QQ(0)], [QQ(1), QQ(0), QQ(0)]] >>> g = [[QQ(1,2)], [QQ(1), QQ(0)], []] >>> dmp_qq_heu_gcd(f, g, 1, QQ) ([[1/1], [2/1, 0/1]], [[1/4], [1/2, 0/1]], [[1/2], []]) """ result = _dmp_ff_trivial_gcd(f, g, u, K0) if result is not None: return result K1 = K0.get_ring() cf, f = dmp_clear_denoms(f, u, K0, K1) cg, g = dmp_clear_denoms(g, u, K0, K1) f = dmp_convert(f, u, K0, K1) g = dmp_convert(g, u, K0, K1) h, cff, cfg = dmp_zz_heu_gcd(f, g, u, K1) h = dmp_convert(h, u, K1, K0) c = dmp_ground_LC(h, u, K0) h = dmp_ground_monic(h, u, K0) cff = dmp_convert(cff, u, K1, K0) cfg = dmp_convert(cfg, u, K1, K0) cff = dmp_mul_ground(cff, K0.quo(c, cf), u, K0) cfg = dmp_mul_ground(cfg, K0.quo(c, cg), u, K0) return h, cff, cfg
def dup_inner_gcd(f, g, K): """ Computes polynomial GCD and cofactors of `f` and `g` in `K[x]`. Returns ``(h, cff, cfg)`` such that ``a = gcd(f, g)``, ``cff = quo(f, h)``, and ``cfg = quo(g, h)``. Examples ======== >>> from sympy.polys.domains import ZZ >>> from sympy.polys.euclidtools import dup_inner_gcd >>> f = ZZ.map([1, 0, -1]) >>> g = ZZ.map([1, -3, 2]) >>> dup_inner_gcd(f, g, ZZ) ([1, -1], [1, 1], [1, -2]) """ if K.has_Field or not K.is_Exact: if K.is_QQ and query('USE_HEU_GCD'): try: return dup_qq_heu_gcd(f, g, K) except HeuristicGCDFailed: pass return dup_ff_prs_gcd(f, g, K) else: if K.is_ZZ and query('USE_HEU_GCD'): try: return dup_zz_heu_gcd(f, g, K) except HeuristicGCDFailed: pass return dup_rr_prs_gcd(f, g, K) @cythonized("u") def _dmp_inner_gcd(f, g, u, K): """Helper function for `dmp_inner_gcd()`. """ if K.has_Field or not K.is_Exact: if K.is_QQ and query('USE_HEU_GCD'): try: return dmp_qq_heu_gcd(f, g, u, K) except HeuristicGCDFailed: pass return dmp_ff_prs_gcd(f, g, u, K) else: if K.is_ZZ and query('USE_HEU_GCD'): try: return dmp_zz_heu_gcd(f, g, u, K) except HeuristicGCDFailed: pass return dmp_rr_prs_gcd(f, g, u, K) @cythonized("u")
[docs]def dmp_inner_gcd(f, g, u, K): """ Computes polynomial GCD and cofactors of `f` and `g` in `K[X]`. Returns ``(h, cff, cfg)`` such that ``a = gcd(f, g)``, ``cff = quo(f, h)``, and ``cfg = quo(g, h)``. Examples ======== >>> from sympy.polys.domains import ZZ >>> from sympy.polys.euclidtools import dmp_inner_gcd >>> f = ZZ.map([[1], [2, 0], [1, 0, 0]]) >>> g = ZZ.map([[1], [1, 0], []]) >>> dmp_inner_gcd(f, g, 1, ZZ) ([[1], [1, 0]], [[1], [1, 0]], [[1], []]) """ if not u: return dup_inner_gcd(f, g, K) J, (f, g) = dmp_multi_deflate((f, g), u, K) h, cff, cfg = _dmp_inner_gcd(f, g, u, K) return (dmp_inflate(h, J, u, K), dmp_inflate(cff, J, u, K), dmp_inflate(cfg, J, u, K))
def dup_gcd(f, g, K): """ Computes polynomial GCD of `f` and `g` in `K[x]`. Examples ======== >>> from sympy.polys.domains import ZZ >>> from sympy.polys.euclidtools import dup_gcd >>> f = ZZ.map([1, 0, -1]) >>> g = ZZ.map([1, -3, 2]) >>> dup_gcd(f, g, ZZ) [1, -1] """ return dup_inner_gcd(f, g, K)[0] @cythonized("u")
[docs]def dmp_gcd(f, g, u, K): """ Computes polynomial GCD of `f` and `g` in `K[X]`. Examples ======== >>> from sympy.polys.domains import ZZ >>> from sympy.polys.euclidtools import dmp_gcd >>> f = ZZ.map([[1], [2, 0], [1, 0, 0]]) >>> g = ZZ.map([[1], [1, 0], []]) >>> dmp_gcd(f, g, 1, ZZ) [[1], [1, 0]] """ return dmp_inner_gcd(f, g, u, K)[0]
def dup_rr_lcm(f, g, K): """ Computes polynomial LCM over a ring in `K[x]`. Examples ======== >>> from sympy.polys.domains import ZZ >>> from sympy.polys.euclidtools import dup_rr_lcm >>> f = ZZ.map([1, 0, -1]) >>> g = ZZ.map([1, -3, 2]) >>> dup_rr_lcm(f, g, ZZ) [1, -2, -1, 2] """ fc, f = dup_primitive(f, K) gc, g = dup_primitive(g, K) c = K.lcm(fc, gc) h = dup_quo(dup_mul(f, g, K), dup_gcd(f, g, K), K) return dup_mul_ground(h, c, K) def dup_ff_lcm(f, g, K): """ Computes polynomial LCM over a field in `K[x]`. Examples ======== >>> from sympy.polys.domains import QQ >>> from sympy.polys.euclidtools import dup_ff_lcm >>> f = [QQ(1,2), QQ(7,4), QQ(3,2)] >>> g = [QQ(1,2), QQ(1), QQ(0)] >>> dup_ff_lcm(f, g, QQ) [1/1, 7/2, 3/1, 0/1] """ h = dup_quo(dup_mul(f, g, K), dup_gcd(f, g, K), K) return dup_monic(h, K) def dup_lcm(f, g, K): """ Computes polynomial LCM of `f` and `g` in `K[x]`. Examples ======== >>> from sympy.polys.domains import ZZ >>> from sympy.polys.euclidtools import dup_lcm >>> f = ZZ.map([1, 0, -1]) >>> g = ZZ.map([1, -3, 2]) >>> dup_lcm(f, g, ZZ) [1, -2, -1, 2] """ if K.has_Field or not K.is_Exact: return dup_ff_lcm(f, g, K) else: return dup_rr_lcm(f, g, K) @cythonized("u") def dmp_rr_lcm(f, g, u, K): """ Computes polynomial LCM over a ring in `K[X]`. Examples ======== >>> from sympy.polys.domains import ZZ >>> from sympy.polys.euclidtools import dmp_rr_lcm >>> f = ZZ.map([[1], [2, 0], [1, 0, 0]]) >>> g = ZZ.map([[1], [1, 0], []]) >>> dmp_rr_lcm(f, g, 1, ZZ) [[1], [2, 0], [1, 0, 0], []] """ fc, f = dmp_ground_primitive(f, u, K) gc, g = dmp_ground_primitive(g, u, K) c = K.lcm(fc, gc) h = dmp_quo(dmp_mul(f, g, u, K), dmp_gcd(f, g, u, K), u, K) return dmp_mul_ground(h, c, u, K) @cythonized("u") def dmp_ff_lcm(f, g, u, K): """ Computes polynomial LCM over a field in `K[X]`. Examples ======== >>> from sympy.polys.domains import QQ >>> from sympy.polys.euclidtools import dmp_ff_lcm >>> f = [[QQ(1,4)], [QQ(1), QQ(0)], [QQ(1), QQ(0), QQ(0)]] >>> g = [[QQ(1,2)], [QQ(1), QQ(0)], []] >>> dmp_ff_lcm(f, g, 1, QQ) [[1/1], [4/1, 0/1], [4/1, 0/1, 0/1], []] """ h = dmp_quo(dmp_mul(f, g, u, K), dmp_gcd(f, g, u, K), u, K) return dmp_ground_monic(h, u, K) @cythonized("u")
[docs]def dmp_lcm(f, g, u, K): """ Computes polynomial LCM of `f` and `g` in `K[X]`. Examples ======== >>> from sympy.polys.domains import ZZ >>> from sympy.polys.euclidtools import dmp_lcm >>> f = ZZ.map([[1], [2, 0], [1, 0, 0]]) >>> g = ZZ.map([[1], [1, 0], []]) >>> dmp_lcm(f, g, 1, ZZ) [[1], [2, 0], [1, 0, 0], []] """ if not u: return dup_lcm(f, g, K) if K.has_Field or not K.is_Exact: return dmp_ff_lcm(f, g, u, K) else: return dmp_rr_lcm(f, g, u, K)
@cythonized("u,v")
[docs]def dmp_content(f, u, K): """ Returns GCD of multivariate coefficients. Examples ======== >>> from sympy.polys.domains import ZZ >>> from sympy.polys.euclidtools import dmp_content >>> f = ZZ.map([[2, 6], [4, 12]]) >>> dmp_content(f, 1, ZZ) [2, 6] """ cont, v = dmp_LC(f, K), u-1 if dmp_zero_p(f, u): return cont for c in f[1:]: cont = dmp_gcd(cont, c, v, K) if dmp_one_p(cont, v, K): break if K.is_negative(dmp_ground_LC(cont, v, K)): return dmp_neg(cont, v, K) else: return cont
@cythonized("u,v")
[docs]def dmp_primitive(f, u, K): """ Returns multivariate content and a primitive polynomial. Examples ======== >>> from sympy.polys.domains import ZZ >>> from sympy.polys.euclidtools import dmp_primitive >>> f = ZZ.map([[2, 6], [4, 12]]) >>> dmp_primitive(f, 1, ZZ) ([2, 6], [[1], [2]]) """ cont, v = dmp_content(f, u, K), u-1 if dmp_zero_p(f, u) or dmp_one_p(cont, v, K): return cont, f else: return cont, [ dmp_quo(c, cont, v, K) for c in f ]
def dup_cancel(f, g, K, include=True): """ Cancel common factors in a rational function `f/g`. Examples ======== >>> from sympy.polys.domains import ZZ >>> from sympy.polys.euclidtools import dup_cancel >>> f = ZZ.map([2, 0, -2]) >>> g = ZZ.map([1, -2, 1]) >>> dup_cancel(f, g, ZZ) ([2, 2], [1, -1]) """ return dmp_cancel(f, g, 0, K, include=include)
[docs]def dmp_cancel(f, g, u, K, include=True): """ Cancel common factors in a rational function `f/g`. Examples ======== >>> from sympy.polys.domains import ZZ >>> from sympy.polys.euclidtools import dmp_cancel >>> f = ZZ.map([[2], [0], [-2]]) >>> g = ZZ.map([[1], [-2], [1]]) >>> dmp_cancel(f, g, 1, ZZ) ([[2], [2]], [[1], [-1]]) """ K0 = None if K.has_Field and K.has_assoc_Ring: K0, K = K, K.get_ring() cq, f = dmp_clear_denoms(f, u, K0, K, convert=True) cp, g = dmp_clear_denoms(g, u, K0, K, convert=True) else: cp, cq = K.one, K.one _, p, q = dmp_inner_gcd(f, g, u, K) if K0 is not None: p = dmp_convert(p, u, K, K0) q = dmp_convert(q, u, K, K0) K = K0 p_neg = K.is_negative(dmp_ground_LC(p, u, K)) q_neg = K.is_negative(dmp_ground_LC(q, u, K)) if p_neg and q_neg: p, q = dmp_neg(p, u, K), dmp_neg(q, u, K) elif p_neg: cp, p = -cp, dmp_neg(p, u, K) elif q_neg: cp, q = -cp, dmp_neg(q, u, K) if not include: return cp, cq, p, q p = dmp_mul_ground(p, cp, u, K) q = dmp_mul_ground(q, cq, u, K) return p, q