(f, g, p, K)
return h[:dq + 1], gf_strip(h[dq + 1:])
Division with remainder in ``GF(p)[x]``.

    Given univariate polynomials ``f`` and ``g`` with coefficients in a
    finite field with ``p`` elements, returns polynomials ``q`` and ``r``
    (quotient and remainder) such that ``f = q*g + r``.

    Consider polynomials ``x**3 + x + 1`` and ``x**2 + x`` in GF(2)::

       >>> from sympy.polys.domains import ZZ
       >>> from sympy.polys.galoistools import gf_div, gf_add_mul

       >>> gf_div(ZZ.map([1, 0, 1, 1]), ZZ.map([1, 1, 0]), 2, ZZ)
       ([1, 1], [1])

    As result we obtained quotient ``x + 1`` and remainder ``1``, thus::

       >>> gf_add_mul(ZZ.map([1]), ZZ.map([1, 1]), ZZ.map([1, 1, 0]), 2, ZZ)
       [1, 0, 1, 1]

    References
    ==========

    1. [Monagan93]_
    2. [Gathen99]_
Division with remainder in ``GF(p)[x]``.
Division with remainder in GF ( p ) [ x ] .
def gf_div(f, g, p, K):
    """
    Division with remainder in ``GF(p)[x]``.

    Given univariate polynomials ``f`` and ``g`` with coefficients in a
    finite field with ``p`` elements, returns polynomials ``q`` and ``r``
    (quotient and remainder) such that ``f = q*g + r``.

    Consider polynomials ``x**3 + x + 1`` and ``x**2 + x`` in GF(2)::

       >>> from sympy.polys.domains import ZZ
       >>> from sympy.polys.galoistools import gf_div, gf_add_mul

       >>> gf_div(ZZ.map([1, 0, 1, 1]), ZZ.map([1, 1, 0]), 2, ZZ)
       ([1, 1], [1])

    As result we obtained quotient ``x + 1`` and remainder ``1``, thus::

       >>> gf_add_mul(ZZ.map([1]), ZZ.map([1, 1]), ZZ.map([1, 1, 0]), 2, ZZ)
       [1, 0, 1, 1]

    References
    ==========

    1. [Monagan93]_
    2. [Gathen99]_

    """
    df = gf_degree(f)
    dg = gf_degree(g)

    if not g:
        raise ZeroDivisionError("polynomial division")
    elif df < dg:
        return [], f

    inv = K.invert(g[0], p)

    h, dq, dr = list(f), df - dg, dg - 1

    for i in xrange(0, df + 1):
        coeff = h[i]

        for j in xrange(max(0, dg - i), min(df - i, dr) + 1):
            coeff -= h[i + j - dg] * g[dg - j]

        if i <= dq:
            coeff *= inv

        h[i] = coeff % p

    return h[:dq + 1], gf_strip(h[dq + 1:])
def gf_div ( f , g , p , K ) : df = gf_degree ( f ) dg = gf_degree ( g ) if not g : raise ZeroDivisionError ( "polynomial division" ) elif df < dg : return [ ] , f inv = K . invert ( g [ 0 ] , p ) h , dq , dr = list ( f ) , df - dg , dg - 1 for i in xrange ( 0 , df + 1 ) : coeff = h [ i ] for j in xrange ( max ( 0 , dg - i ) , min ( df - i , dr ) + 1 ) : coeff -= h [ i + j - dg ] * g [ dg - j ] if i <= dq : coeff *= inv h [ i ] = coeff % p return h [ : dq + 1 ] , gf_strip ( h [ dq + 1 : ] )