def dup_extract(f, g, K):
    """
    Extract common content from a pair of polynomials in ``K[x]``.

    Examples
    ========

    >>> from sympy.polys.domains import ZZ
    >>> from sympy.polys.densetools import dup_extract

    >>> f = ZZ.map([6, 12, 18])
    >>> g = ZZ.map([4, 8, 12])

    >>> dup_extract(f, g, ZZ)
    (2, [3, 6, 9], [2, 4, 6])

    """
    fc = dup_content(f, K)
    gc = dup_content(g, K)

    gcd = K.gcd(fc, gc)

    if not K.is_one(gcd):
        f = dup_quo_ground(f, gcd, K)
        g = dup_quo_ground(g, gcd, K)

    return gcd, f, g

def dup_primitive(f, K):
    """
    Compute content and the primitive form of ``f`` in ``K[x]``.

    Examples
    ========

    >>> from sympy.polys import ring, ZZ, QQ

    >>> R, x = ring("x", ZZ)
    >>> f = 6*x**2 + 8*x + 12

    >>> R.dup_primitive(f)
    (2, 3*x**2 + 4*x + 6)

    >>> R, x = ring("x", QQ)
    >>> f = 6*x**2 + 8*x + 12

    >>> R.dup_primitive(f)
    (2, 3*x**2 + 4*x + 6)

    """
    if not f:
        return K.zero, f

    cont = dup_content(f, K)

    if K.is_one(cont):
        return cont, f
    else:
        return cont, dup_quo_ground(f, cont, K)

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

def dup_primitive(f, K):
    """
    Compute content and the primitive form of ``f`` in ``K[x]``.

    Examples
    ========

    >>> from sympy.polys.domains import ZZ, QQ
    >>> from sympy.polys.densetools import dup_primitive

    >>> f = ZZ.map([6, 8, 12])
    >>> g = QQ.map([6, 8, 12])

    >>> dup_primitive(f, ZZ)
    (2, [3, 4, 6])
    >>> dup_primitive(g, QQ)
    (2/1, [3/1, 4/1, 6/1])

    """
    if not f:
        return K.zero, f

    cont = dup_content(f, K)

    if K.is_one(cont):
        return cont, f
    else:
        return cont, dup_quo_ground(f, cont, K)

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 import ring, QQ
    >>> R, x = ring("x", QQ)

    >>> f = x**4 - 2*x**3 - 6*x**2 + 12*x + 15
    >>> g = x**3 + x**2 - 4*x - 4

    >>> R.dup_half_gcdex(f, g)
    (-1/5*x + 3/5, x + 1)

    """
    if not K.has_Field:
        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

def dup_factor_list(f, K0):
    """Factor univariate polynomials into irreducibles in `K[x]`. """
    j, f = dup_terms_gcd(f, K0)
    cont, f = dup_primitive(f, K0)

    if K0.is_FiniteField:
        coeff, factors = dup_gf_factor(f, K0)
    elif K0.is_Algebraic:
        coeff, factors = dup_ext_factor(f, K0)
    else:
        if not K0.is_Exact:
            K0_inexact, K0 = K0, K0.get_exact()
            f = dup_convert(f, K0_inexact, K0)
        else:
            K0_inexact = None

        if K0.is_Field:
            K = K0.get_ring()

            denom, f = dup_clear_denoms(f, K0, K)
            f = dup_convert(f, K0, K)
        else:
            K = K0

        if K.is_ZZ:
            coeff, factors = dup_zz_factor(f, K)
        elif K.is_Poly:
            f, u = dmp_inject(f, 0, K)

            coeff, factors = dmp_factor_list(f, u, K.dom)

            for i, (f, k) in enumerate(factors):
                factors[i] = (dmp_eject(f, u, K), k)

            coeff = K.convert(coeff, K.dom)
        else:  # pragma: no cover
            raise DomainError('factorization not supported over %s' % K0)

        if K0.is_Field:
            for i, (f, k) in enumerate(factors):
                factors[i] = (dup_convert(f, K, K0), k)

            coeff = K0.convert(coeff, K)
            coeff = K0.quo(coeff, denom)

            if K0_inexact:
                for i, (f, k) in enumerate(factors):
                    max_norm = dup_max_norm(f, K0)
                    f = dup_quo_ground(f, max_norm, K0)
                    f = dup_convert(f, K0, K0_inexact)
                    factors[i] = (f, k)
                    coeff = K0.mul(coeff, K0.pow(max_norm, k))

                coeff = K0_inexact.convert(coeff, K0)
                K0 = K0_inexact

    if j:
        factors.insert(0, ([K0.one, K0.zero], j))

    return coeff*cont, _sort_factors(factors)

def dup_monic(f, K):
    """
    Divides all coefficients by ``LC(f)`` in ``K[x]``.

    **Examples**

    >>> from sympy.polys.domains import ZZ, QQ
    >>> from sympy.polys.densetools import dup_monic

    >>> dup_monic([ZZ(3), ZZ(6), ZZ(9)], ZZ)
    [1, 2, 3]

    >>> dup_monic([QQ(3), QQ(4), QQ(2)], QQ)
    [1/1, 4/3, 2/3]

    """
    if not f:
        return f

    lc = dup_LC(f, K)

    if K.is_one(lc):
        return f
    else:
        return dup_quo_ground(f, lc, K)

def dup_extract(f, g, K):
    """
    Extract common content from a pair of polynomials in ``K[x]``.

    Examples
    ========

    >>> from sympy.polys import ring, ZZ
    >>> R, x = ring("x", ZZ)

    >>> R.dup_extract(6*x**2 + 12*x + 18, 4*x**2 + 8*x + 12)
    (2, 3*x**2 + 6*x + 9, 2*x**2 + 4*x + 6)

    """
    fc = dup_content(f, K)
    gc = dup_content(g, K)

    gcd = K.gcd(fc, gc)

    if not K.is_one(gcd):
        f = dup_quo_ground(f, gcd, K)
        g = dup_quo_ground(g, gcd, K)

    return gcd, f, g

def test_dup_quo_ground():
    raises(ZeroDivisionError, lambda: dup_quo_ground(dup_normal([1,2,3], ZZ), ZZ(0), ZZ))

    f = dup_normal([], ZZ)

    assert dup_quo_ground(f, ZZ(3), ZZ) == dup_normal([], ZZ)

    f = dup_normal([6,2,8], ZZ)

    assert dup_quo_ground(f, ZZ(1), ZZ) == f
    assert dup_quo_ground(f, ZZ(2), ZZ) == dup_normal([3,1,4], ZZ)

    assert dup_quo_ground(f, ZZ(3), ZZ) == dup_normal([2,0,2], ZZ)

    f = dup_normal([6,2,8], QQ)

    assert dup_quo_ground(f, QQ(1), QQ) == f
    assert dup_quo_ground(f, QQ(2), QQ) == [QQ(3),QQ(1),QQ(4)]
    assert dup_quo_ground(f, QQ(7), QQ) == [QQ(6,7),QQ(2,7),QQ(8,7)]

def test_dup_quo_ground():
    raises(ZeroDivisionError, 'dup_quo_ground(dup_normal([1,2,3], ZZ), ZZ(0), ZZ)')
    raises(ExactQuotientFailed, 'dup_quo_ground(dup_normal([1,2,3], ZZ), ZZ(3), ZZ)')

    f = dup_normal([], ZZ)

    assert dup_quo_ground(f, ZZ(3), ZZ) == dup_normal([], ZZ)

    f = dup_normal([6,2,8], ZZ)

    assert dup_quo_ground(f, ZZ(1), ZZ) == f
    assert dup_quo_ground(f, ZZ(2), ZZ) == dup_normal([3,1,4], ZZ)

    f = dup_normal([6,2,8], QQ)

    assert dup_quo_ground(f, QQ(1), QQ) == f
    assert dup_quo_ground(f, QQ(2), QQ) == [QQ(3),QQ(1),QQ(4)]
    assert dup_quo_ground(f, QQ(7), QQ) == [QQ(6,7),QQ(2,7),QQ(8,7)]

def dup_inner_subresultants(f, g, K):
    """
    Subresultant PRS algorithm in `K[x]`.

    Computes the subresultant polynomial remainder sequence (PRS)
    and the non-zero scalar subresultants of `f` and `g`.
    By [1] Thm. 3, these are the constants '-c' (- to optimize
    computation of sign).
    The first subdeterminant is set to 1 by convention to match
    the polynomial and the scalar subdeterminants.
    If 'deg(f) < deg(g)', the subresultants of '(g,f)' are computed.

    Examples
    ========

    >>> from sympy.polys import ring, ZZ
    >>> R, x = ring("x", ZZ)

    >>> R.dup_inner_subresultants(x**2 + 1, x**2 - 1)
    ([x**2 + 1, x**2 - 1, -2], [1, 1, 4])

    References
    ==========

    [1] W.S. Brown, The Subresultant PRS Algorithm.
    ACM Transaction of Mathematical Software 4 (1978) 237-249

    """
    n = dup_degree(f)
    m = dup_degree(g)

    if n < m:
        f, g = g, f
        n, m = m, n

    if not f:
        return [], []

    if not g:
        return [f], [K.one]

    R = [f, g]
    d = n - m

    b = (-K.one)**(d + 1)

    h = dup_prem(f, g, K)
    h = dup_mul_ground(h, b, K)

    lc = dup_LC(g, K)
    c = lc**d

    # Conventional first scalar subdeterminant is 1
    S = [K.one, c]
    c = -c

    while h:
        k = dup_degree(h)
        R.append(h)

        f, g, m, d = g, h, k, m - k

        b = -lc * c**d

        h = dup_prem(f, g, K)
        h = dup_quo_ground(h, b, K)

        lc = dup_LC(g, K)

        if d > 1:        # abnormal case
            q = c**(d - 1)
            c = K.quo((-lc)**d, q)
        else:
            c = -lc

        S.append(-c)

    return R, S

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