'Division algorithm, see [CLO] p64.
fv array of polynomials
return qv, r such that
self = sum(fv[i]*qv[i]) + r
All polynomials are required not to be Laurent polynomials.
Examples
>>> from sympy.polys.rings import ring
>>> from sympy.polys.domains import ZZ
>>> _, x, y = ring(\'x, y\', ZZ)
>>> f = x**3
>>> f0 = x - y**2
>>> f1 = x - y
>>> qv, r = f.div((f0, f1))
>>> qv[0]
x**2 + x*y**2 + y**4
>>> qv[1]
0
>>> r
y**6'

def div(self, fv):

    ring = self.ring
    domain = ring.domain
    ret_single = False
    if isinstance(fv, PolyElement):
        ret_single = True
        fv = [fv]
    if any(((not f) for f in fv)):
        raise ZeroDivisionError('polynomial    division')
    if (not self):
        if ret_single:
            return (ring.zero, ring.zero)
        else:
            return ([], ring.zero)
    for f in fv:
        if (f.ring != ring):
            raise ValueError('self    and    f    must    have    the    same    ring')
    s = len(fv)
    qv = [ring.zero for i in range(s)]
    p = self.copy()
    r = ring.zero
    term_div = self._term_div()
    expvs = [fx.leading_expv() for fx in fv]
    while p:
        i = 0
        divoccurred = 0
        while ((i < s) and (divoccurred == 0)):
            expv = p.leading_expv()
            term = term_div((expv, p[expv]), (expvs[i], fv[i][expvs[i]]))
            if (term is not None):
                (expv1, c) = term
                qv[i] = qv[i]._iadd_monom((expv1, c))
                p = p._iadd_poly_monom(fv[i], (expv1, (- c)))
                divoccurred = 1
            else:
                i += 1
        if (not divoccurred):
            expv = p.leading_expv()
            r = r._iadd_monom((expv, p[expv]))
            del p[expv]
    if (expv == ring.zero_monom):
        r += p
    if ret_single:
        if (not qv):
            return (ring.zero, r)
        else:
            return (qv[0], r)
    else:
        return (qv, r)
