def _eval_rewrite_as_Add(self, *args, **kwargs):
        """return Eq(L, R) as L - R. To control the evaluation of
        the result set pass `evaluate=True` to give L - R;
        if `evaluate=None` then terms in L and R will not cancel
        but they will be listed in canonical order; otherwise
        non-canonical args will be returned.

        Examples
        ========

        >>> from sympy import Eq, Add
        >>> from sympy.abc import b, x
        >>> eq = Eq(x + b, x - b)
        >>> eq.rewrite(Add)
        2*b
        >>> eq.rewrite(Add, evaluate=None).args
        (b, b, x, -x)
        >>> eq.rewrite(Add, evaluate=False).args
        (b, x, b, -x)
        """
        L, R = args
        evaluate = kwargs.get('evaluate', True)
        if evaluate:
            # allow cancellation of args
            return L - R
        args = Add.make_args(L) + Add.make_args(-R)
        if evaluate is None:
            # no cancellation, but canonical
            return _unevaluated_Add(*args)
        # no cancellation, not canonical
        return Add._from_args(args)

def _gcd_terms(terms, isprimitive=False, fraction=True):
    """Helper function for :func:`gcd_terms`.

    If ``isprimitive`` is True then the call to primitive
    for an Add will be skipped. This is useful when the
    content has already been extrated.

    If ``fraction`` is True then the expression will appear over a common
    denominator, the lcm of all term denominators.
    """

    if isinstance(terms, Basic) and not isinstance(terms, Tuple):
        terms = Add.make_args(terms)

    terms = map(Term, [t for t in terms if t])

    # there is some simplification that may happen if we leave this
    # here rather than duplicate it before the mapping of Term onto
    # the terms
    if len(terms) == 0:
        return S.Zero, S.Zero, S.One


    if len(terms) == 1:
        cont = terms[0].coeff
        numer = terms[0].numer.as_expr()
        denom = terms[0].denom.as_expr()

    else:
        cont = terms[0]
        for term in terms[1:]:
            cont = cont.gcd(term)

        for i, term in enumerate(terms):
            terms[i] = term.quo(cont)

        if fraction:
            denom = terms[0].denom

            for term in terms[1:]:
                denom = denom.lcm(term.denom)

            numers = []
            for term in terms:
                numer = term.numer.mul(denom.quo(term.denom))
                numers.append(term.coeff*numer.as_expr())
        else:
            numers = [t.as_expr() for t in terms]
            denom = Term(S(1)).numer

        cont = cont.as_expr()
        numer = Add(*numers)
        denom = denom.as_expr()

    if not isprimitive and numer.is_Add:
        _cont, numer = numer.primitive()
        cont *= _cont

    return cont, numer, denom

def factor_terms(expr, radical=False):
    """Remove common factors from terms in all arguments without
    changing the underlying structure of the expr. No expansion or
    simplification (and no processing of non-commutatives) is performed.

    If radical=True then a radical common to all terms will be factored
    out of any Add sub-expressions of the expr.

    Examples
    ========

    >>> from sympy import factor_terms, Symbol
    >>> from sympy.abc import x, y
    >>> factor_terms(x + x*(2 + 4*y)**3)
    x*(8*(2*y + 1)**3 + 1)
    >>> A = Symbol('A', commutative=False)
    >>> factor_terms(x*A + x*A + x*y*A)
    x*(y*A + 2*A)

    """

    expr = sympify(expr)
    is_iterable = iterable(expr)

    if not isinstance(expr, Basic) or expr.is_Atom:
        if is_iterable:
            return type(expr)([factor_terms(i, radical=radical) for i in expr])
        return expr

    if expr.is_Function or is_iterable or not hasattr(expr, 'args_cnc'):
        args = expr.args
        newargs = tuple([factor_terms(i, radical=radical) for i in args])
        if newargs == args:
            return expr
        return expr.func(*newargs)

    cont, p = expr.as_content_primitive(radical=radical)
    list_args, nc = zip(*[ai.args_cnc() for ai in Add.make_args(p)])
    list_args = list(list_args)
    nc = [((Dummy(), Mul._from_args(i)) if i else None) for i in nc]
    ncreps = dict([i for i in nc if i is not None])
    for i, a in enumerate(list_args):
        if nc[i] is not None:
            a.append(nc[i][0])
        a = Mul._from_args(a) # gcd_terms will fix up ordering
        list_args[i] = gcd_terms(a, isprimitive=True)
        # cancel terms that may not have cancelled
    p = Add._from_args(list_args) # gcd_terms will fix up ordering
    p = gcd_terms(p, isprimitive=True).xreplace(ncreps)
    return _keep_coeff(cont, p)

    def _eval_imageset(self, f):
        expr = f.expr
        if not isinstance(expr, Expr):
            return

        if len(f.variables) > 1:
            return

        n = f.variables[0]

        # f(x) + c and f(-x) + c cover the same integers
        # so choose the form that has the fewest negatives
        c = f(0)
        fx = f(n) - c
        f_x = f(-n) - c
        neg_count = lambda e: sum(_coeff_isneg(_) for _ in Add.make_args(e))
        if neg_count(f_x) < neg_count(fx):
            expr = f_x + c

        a = Wild('a', exclude=[n])
        b = Wild('b', exclude=[n])
        match = expr.match(a*n + b)
        if match and match[a]:
            # canonical shift
            expr = match[a]*n + match[b] % match[a]

        if expr != f.expr:
            return ImageSet(Lambda(n, expr), S.Integers)

    def _eval_simplify(self, ratio=1.7, measure=None):
        from sympy.simplify.simplify import factor_sum, sum_combine
        from sympy.core.function import expand
        from sympy.core.mul import Mul

        # split the function into adds
        terms = Add.make_args(expand(self.function))
        s_t = [] # Sum Terms
        o_t = [] # Other Terms

        for term in terms:
            if term.has(Sum):
                # if there is an embedded sum here
                # it is of the form x * (Sum(whatever))
                # hence we make a Mul out of it, and simplify all interior sum terms
                subterms = Mul.make_args(expand(term))
                out_terms = []
                for subterm in subterms:
                    # go through each term
                    if isinstance(subterm, Sum):
                        # if it's a sum, simplify it
                        out_terms.append(subterm._eval_simplify())
                    else:
                        # otherwise, add it as is
                        out_terms.append(subterm)

                # turn it back into a Mul
                s_t.append(Mul(*out_terms))
            else:
                o_t.append(term)

        # next try to combine any interior sums for further simplification
        result = Add(sum_combine(s_t), *o_t)

        return factor_sum(result, limits=self.limits)

def _transform_explike_DE(DE, g, x, order, syms):
    """Converts DE with free parameters into DE with constant coefficients."""
    from sympy.solvers.solveset import linsolve

    eq = []
    highest_coeff = DE.coeff(Derivative(g(x), x, order))
    for i in range(order):
        coeff = DE.coeff(Derivative(g(x), x, i))
        coeff = (coeff / highest_coeff).expand().collect(x)
        for t in Add.make_args(coeff):
            eq.append(t)
    temp = []
    for e in eq:
        if e.has(x):
            break
        elif e.has(Symbol):
            temp.append(e)
    else:
        eq = temp
    if eq:
        sol = dict(zip(syms, (i for s in linsolve(eq, list(syms)) for i in s)))
        if sol:
            DE = DE.subs(sol)
            DE = DE.factor().as_coeff_mul(Derivative)[1][0]
            DE = DE.collect(Derivative(g(x)))
    return DE

def _peeloff_pi(arg):
    """
    Split ARG into two parts, a "rest" and a multiple of pi/2.
    This assumes ARG to be an Add.
    The multiple of pi returned in the second position is always a Rational.

    Examples:
    >>> from sympy.functions.elementary.trigonometric import _peeloff_pi as peel
    >>> from sympy import pi
    >>> from sympy.abc import x, y
    >>> peel(x + pi/2)
    (x, pi/2)
    >>> peel(x + 2*pi/3 + pi*y)
    (x + pi*y + pi/6, pi/2)
    """
    for a in Add.make_args(arg):
        if a is S.Pi:
            K = S.One
            break
        elif a.is_Mul:
            K, p = a.as_two_terms()
            if p is S.Pi and K.is_Rational:
                break
    else:
        return arg, S.Zero

    m1 = (K % S.Half) * S.Pi
    m2 = K*S.Pi - m1
    return arg - m2, m2

def _peeloff_ipi(arg):
    """
    Split ARG into two parts, a "rest" and a multiple of I*pi/2.
    This assumes ARG to be an Add.
    The multiple of I*pi returned in the second position is always a Rational.

    Examples
    ========

    >>> from sympy.functions.elementary.hyperbolic import _peeloff_ipi as peel
    >>> from sympy import pi, I
    >>> from sympy.abc import x, y
    >>> peel(x + I*pi/2)
    (x, I*pi/2)
    >>> peel(x + I*2*pi/3 + I*pi*y)
    (x + I*pi*y + I*pi/6, I*pi/2)
    """
    for a in Add.make_args(arg):
        if a == S.Pi*S.ImaginaryUnit:
            K = S.One
            break
        elif a.is_Mul:
            K, p = a.as_two_terms()
            if p == S.Pi*S.ImaginaryUnit and K.is_Rational:
                break
    else:
        return arg, S.Zero

    m1 = (K % S.Half)*S.Pi*S.ImaginaryUnit
    m2 = K*S.Pi*S.ImaginaryUnit - m1
    return arg - m2, m2

    def do(expr):
        from sympy.concrete.summations import Sum
        from sympy.simplify.simplify import factor_sum
        is_iterable = iterable(expr)

        if not isinstance(expr, Basic) or expr.is_Atom:
            if is_iterable:
                return type(expr)([do(i) for i in expr])
            return expr

        if expr.is_Pow or expr.is_Function or \
                is_iterable or not hasattr(expr, 'args_cnc'):
            args = expr.args
            newargs = tuple([do(i) for i in args])
            if newargs == args:
                return expr
            return expr.func(*newargs)

        if isinstance(expr, Sum):
            return factor_sum(expr, radical=radical, clear=clear, fraction=fraction, sign=sign)

        cont, p = expr.as_content_primitive(radical=radical, clear=clear)
        if p.is_Add:
            list_args = [do(a) for a in Add.make_args(p)]
            # get a common negative (if there) which gcd_terms does not remove
            if all(a.as_coeff_Mul()[0].extract_multiplicatively(-1) is not None
                   for a in list_args):
                cont = -cont
                list_args = [-a for a in list_args]
            # watch out for exp(-(x+2)) which gcd_terms will change to exp(-x-2)
            special = {}
            for i, a in enumerate(list_args):
                b, e = a.as_base_exp()
                if e.is_Mul and e != Mul(*e.args):
                    list_args[i] = Dummy()
                    special[list_args[i]] = a
            # rebuild p not worrying about the order which gcd_terms will fix
            p = Add._from_args(list_args)
            p = gcd_terms(p,
                isprimitive=True,
                clear=clear,
                fraction=fraction).xreplace(special)
        elif p.args:
            p = p.func(
                *[do(a) for a in p.args])
        rv = _keep_coeff(cont, p, clear=clear, sign=sign)
        return rv

def factor_terms(expr):
    """Remove common factors from terms in all arguments without
    changing the underlying structure of the expr. No expansion or
    simplification (and no processing of non-commutative) is performed.

    **Examples**

    >>> from sympy import factor_terms, Symbol
    >>> from sympy.abc import x, y
    >>> factor_terms(x + x*(2 + 4*y)**3)
    x*(8*(2*y + 1)**3 + 1)
    >>> A = Symbol('A', commutative=False)
    >>> factor_terms(x*A + x*A + x*y*A)
    x*(y*A + 2*A)

    """

    expr = sympify(expr)

    if iterable(expr):
        return type(expr)([factor_terms(i) for i in expr])

    if not isinstance(expr, Basic) or expr.is_Atom:
        return expr

    if expr.is_Function:
        return expr.func(*[factor_terms(i) for i in expr.args])

    cont, p = expr.as_content_primitive()
    list_args, nc = zip(*[ai.args_cnc(clist=True) for ai in Add.make_args(p)])
    list_args = list(list_args)
    nc = [((Dummy(), Mul._from_args(i)) if i else None) for i in nc]
    ncreps = dict([i for i in nc if i is not None])
    for i, a in enumerate(list_args):
        if nc[i] is not None:
           a.append(nc[i][0])
        a = Mul._from_args(a) # gcd_terms will fix up ordering
        list_args[i] = gcd_terms(a, isprimitive=True)
        # cancel terms that may not have cancelled
    p = Add._from_args(list_args) # gcd_terms will fix up ordering
    p = gcd_terms(p, isprimitive=True).subs(ncreps) # exact subs could be used here
    return _keep_coeff(cont, p)

def _gcd_terms(terms, isprimitive=False):
    """Helper function for :func:`gcd_terms`. If `isprimitive` is True then the
    call to primitive for an Add will be skipped. This is useful when the
    content has already been extrated."""
    if isinstance(terms, Basic) and not isinstance(terms, Tuple):
        terms = Add.make_args(terms)

    if len(terms) <= 1:
        if not terms:
            return S.Zero, S.Zero, S.One
        else:
            return terms[0], S.One, S.One

    terms = map(Term, terms)
    cont = terms[0]

    for term in terms[1:]:
        cont = cont.gcd(term)

    for i, term in enumerate(terms):
        terms[i] = term.quo(cont)

    denom = terms[0].denom

    for term in terms[1:]:
        denom = denom.lcm(term.denom)

    numers = []

    for term in terms:
        numer = term.numer.mul(denom.quo(term.denom))
        numers.append(term.coeff*numer.as_expr())

    cont = cont.as_expr()
    numer = Add(*numers)
    denom = denom.as_expr()
    if not isprimitive and numer.is_Add:
        _cont, numer = numer.primitive()
        cont *= _cont

    return cont, numer, denom

def _solve_hyper_RE(f, x, RE, g, k):
    """See docstring of :func:`rsolve_hypergeometric` for details."""
    terms = Add.make_args(RE)

    if len(terms) == 2:
        gs = list(RE.atoms(Function))
        P, Q = map(RE.coeff, gs)
        m = gs[1].args[0] - gs[0].args[0]
        if m < 0:
            P, Q = Q, P
            m = abs(m)
        return rsolve_hypergeometric(f, x, P, Q, k, m)

def _solve_explike_DE(f, x, DE, g, k):
    """Solves DE with constant coefficients."""
    from sympy.solvers import rsolve

    for t in Add.make_args(DE):
        coeff, d = t.as_independent(g)
        if coeff.free_symbols:
            return

    RE = exp_re(DE, g, k)

    init = {}
    for i in range(len(Add.make_args(RE))):
        if i:
            f = f.diff(x)
        init[g(k).subs(k, i)] = f.limit(x, 0)

    sol = rsolve(RE, g(k), init)

    if sol:
        return (sol / factorial(k), S.Zero, S.Zero)

def _gcd_terms(terms):
    """Helper function for :func:`gcd_terms`. """
    if isinstance(terms, Basic):
        terms = Add.make_args(terms)

    if len(terms) <= 1:
        if not terms:
            return S.Zero, S.Zero, S.One
        else:
            return terms[0], S.One, S.One

    terms = map(Term, terms)
    cont = terms[0]

    for term in terms[1:]:
        cont = cont.gcd(term)

    for i, term in enumerate(terms):
        terms[i] = term.quo(cont)

    denom = terms[0].denom

    for term in terms[1:]:
        denom = denom.lcm(term.denom)

    numers = []

    for term in terms:
        numer = term.numer.mul(denom.quo(term.denom))
        numers.append(term.coeff*numer.as_expr())

    cont = cont.as_expr()
    numer = Add(*numers)
    denom = denom.as_expr()

    if numer.is_Add:
        _cont, numer = numer.primitive()
        cont *= _cont

    return cont, numer, denom

def hyper_re(DE, r, k):
    """Converts a DE into a RE.

    Performs the substitution:

    .. math::
        x^l f^j(x) \\to (k + 1 - l)_j . a_{k + j - l}

    Normalises the terms so that lowest order of a term is always r(k).

    Examples
    ========

    >>> from sympy import Function, Derivative
    >>> from sympy.series.formal import hyper_re
    >>> from sympy.abc import x, k
    >>> f, r = Function('f'), Function('r')

    >>> hyper_re(-f(x) + Derivative(f(x)), r, k)
    (k + 1)*r(k + 1) - r(k)
    >>> hyper_re(-x*f(x) + Derivative(f(x), x, x), r, k)
    (k + 2)*(k + 3)*r(k + 3) - r(k)

    See Also
    ========

    sympy.series.formal.exp_re
    """
    RE = S.Zero

    g = DE.atoms(Function).pop()
    x = g.atoms(Symbol).pop()

    mini = None
    for t in Add.make_args(DE.expand()):
        coeff, d = t.as_independent(g)
        c, v = coeff.as_independent(x)
        l = v.as_coeff_exponent(x)[1]
        if isinstance(d, Derivative):
            j = len(d.args[1:])
        else:
            j = 0
        RE += c * rf(k + 1 - l, j) * r(k + j - l)
        if mini is None or j - l < mini:
            mini = j - l

    RE = RE.subs(k, k - mini)

    m = Wild('m')
    return RE.collect(r(k + m))

    def calc_part(expr, nexpr):
        from sympy.core.add import Add
        nint = int(to_int(nexpr, rnd))
        n, c, p, b = nexpr
        is_int = (p == 0)
        if not is_int:
            # if there are subs and they all contain integer re/im parts
            # then we can (hopefully) safely substitute them into the
            # expression
            s = options.get('subs', False)
            if s:
                doit = True
                from sympy.core.compatibility import as_int
                for v in s.values():
                    try:
                        as_int(v)
                    except ValueError:
                        try:
                            [as_int(i) for i in v.as_real_imag()]
                            continue
                        except (ValueError, AttributeError):
                            doit = False
                            break
                if doit:
                    expr = expr.subs(s)

            expr = Add(expr, -nint, evaluate=False)
            x, _, x_acc, _ = evalf(expr, 10, options)
            try:
                check_target(expr, (x, None, x_acc, None), 3)
            except PrecisionExhausted:
                if not expr.equals(0):
                    raise PrecisionExhausted
                x = fzero
            nint += int(no*(mpf_cmp(x or fzero, fzero) == no))
        nint = from_int(nint)
        return nint, fastlog(nint) + 10

def _solve_simple(f, x, DE, g, k):
    """Converts DE into RE and solves using :func:`rsolve`."""
    from sympy.solvers import rsolve

    RE = hyper_re(DE, g, k)

    init = {}
    for i in range(len(Add.make_args(RE))):
        if i:
            f = f.diff(x)
        init[g(k).subs(k, i)] = f.limit(x, 0) / factorial(i)

    sol = rsolve(RE, g(k), init)

    if sol:
        return (sol, S.Zero, S.Zero)

def exp_re(DE, r, k):
    """Converts a DE with constant coefficients (explike) into a RE.

    Performs the substitution:

    .. math::
        f^j(x) \\to r(k + j)

    Normalises the terms so that lowest order of a term is always r(k).

    Examples
    ========

    >>> from sympy import Function, Derivative
    >>> from sympy.series.formal import exp_re
    >>> from sympy.abc import x, k
    >>> f, r = Function('f'), Function('r')

    >>> exp_re(-f(x) + Derivative(f(x)), r, k)
    -r(k) + r(k + 1)
    >>> exp_re(Derivative(f(x), x) + Derivative(f(x), x, x), r, k)
    r(k) + r(k + 1)

    See Also
    ========

    sympy.series.formal.hyper_re
    """
    RE = S.Zero

    g = DE.atoms(Function).pop()

    mini = None
    for t in Add.make_args(DE):
        coeff, d = t.as_independent(g)
        if isinstance(d, Derivative):
            j = len(d.args) - 1
        else:
            j = 0
        if mini is None or j < mini:
            mini = j
        RE += coeff * r(k + j)
    if mini:
        RE = RE.subs(k, k - mini)
    return RE

def _transform_explike_DE(DE, g, x, order, syms):
    """Converts DE with free parameters into DE with constant coefficients."""
    from sympy.solvers import solve

    eq = []
    highest_coeff = DE.coeff(Derivative(g(x), x, order))
    for i in range(order):
        coeff = DE.coeff(Derivative(g(x), x, i))
        coeff = (coeff / highest_coeff).expand().collect(x)
        for t in Add.make_args(coeff):
            if t.has(x):
                eq.append(t)
    sol = solve(eq, syms, dict=True)
    if sol:
        DE = DE.subs(sol[0])
        DE = DE.factor().as_coeff_mul(Derivative)[1][0]
        DE = DE.collect(Derivative(g(x)))
    return DE

    def integrate(field=None):
        irreducibles = set()

        for poly in reducibles:
            for z in poly.atoms(Symbol):
                if z in V:
                    break
            else:
                continue

            irreducibles |= set(root_factors(poly, z, filter=field))

        log_coeffs, log_part = [], []
        B = _symbols('B', len(irreducibles))

        for i, poly in enumerate(irreducibles):
            if poly.has(*V):
                log_coeffs.append(B[i])
                log_part.append(log_coeffs[-1] * log(poly))

        coeffs = poly_coeffs + log_coeffs

        candidate = poly_part/poly_denom + Add(*log_part)

        h = F - derivation(candidate) / denom

        numer = h.as_numer_denom()[0].expand()

        equations = {}

        for term in Add.make_args(numer):
            coeff, dependent = term.as_independent(*V)

            if dependent in equations:
                equations[dependent] += coeff
            else:
                equations[dependent] = coeff

        solution = solve(equations.values(), *coeffs)

        if solution is not None:
            return (solution, candidate, coeffs)
        else:
            return None