How to use the sympy.Add function in sympy

# Leads to 2D aliases
        eqn = Eq(u.forward, ((u[t, x, y, z] + u[t, x, y+1, z+1])*3*f +
                             (u[t, x, y+2, z+2] + u[t, x, y+3, z+3])*3*f + 1))
        op0 = Operator(eqn, opt=('noop', {'openmp': True}))
        op1 = Operator(eqn, opt=('advanced', {'openmp': True, 'cire-mincost-sops': 1}))

        y0_blk_size = op1.parameters[2]
        z_size = op1.parameters[3]

        arrays = [i for i in FindSymbols().visit(op1._func_table['bf0'].root)
                  if i.is_Array and i._mem_local]
        assert len(arrays) == 1
        a = arrays[0]
        assert len(a.dimensions) == 2
        assert a.halo == ((1, 1), (1, 1))
        assert Add(*a.symbolic_shape[0].args) == y0_blk_size + 2
        assert Add(*a.symbolic_shape[1].args) == z_size + 2

        # Check numerical output
        op0(time_M=1)
        exp = np.copy(u.data[:])
        u.data_with_halo[:] = 0.5
        op1(time_M=1)
        assert np.all(u.data == exp)

def testComplex(self):
        self.compare(
            mathics.Complex(mathics.Real('1.0'), mathics.Real('1.0')),
            sympy.Add(sympy.Float('1.0'), sympy.Float('1.0') * sympy.I))

        self.compare(
            mathics.Complex(mathics.Integer(0), mathics.Integer(1)),
            sympy.I)

        self.compare(
            mathics.Complex(mathics.Integer(-1), mathics.Integer(1)),
            sympy.Integer(-1) + sympy.I)

if complex_pair:
            a, b = e.as_real_imag()
            if b:
                return "(%s;%s)" % (encode_item(*expr_form(a)),
                                    encode_item(*expr_form(b))), "s"
        a, b = e.as_coeff_mul()
        if not b:
            return number_form(a)
        if abs(a.p) != 1 and abs(sp.Mul(*b).as_numer_denom()[0]) != 1:
            return op("*", a, sp.Mul(*b))
        if a == -1 and b:
            return func("neg", sp.Mul(*b))
        a, b = e.as_coeff_add()
        if a != 0 and b:
            if b[0].as_coeff_mul()[0] < 0:
                return op("-", a, -sp.Add(*b))
            else:
                return op("+", a, sp.Add(*b))
        a, b = e.as_numer_denom()
        if b != 1:
            for r in a.as_coeff_mul()[1]:
                r, p = r.as_base_exp()
                if p.is_Rational and p.p == 1 and p.q > 1:
                    r = r.as_coeff_mul()[0]
                    if r.q == 1 and b.as_coeff_mul()[0] % r.p == 0:
                        a, b = a / r.p**p, b / r.p**p
            if b.as_coeff_mul()[0] != 1:
                if a.as_coeff_mul()[0] == -1 or all(
                        x.as_coeff_mul()[0]<0 for x in a.as_coeff_add()[1]):
                    a, b = -a, -b
            return op("%", a, b)
        if e.func == sp.Add:

def expr_to_json(expr):
    """
    Converts a Sympy expression to a json-compatible tree-structure.
    """
    if isinstance(expr, sympy.Mul):
        return {"type": "Mul", "args": [expr_to_json(arg) for arg in expr.args]}
    elif isinstance(expr, sympy.Add):
        return {"type": "Add", "args": [expr_to_json(arg) for arg in expr.args]}
    elif isinstance(expr, sympy.Symbol):
        return {"type": "Symbol", "name": expr.name}
    elif isinstance(expr, sympy.Pow):
        return {"type": "Pow", "args": [expr_to_json(arg) for arg in expr.args]}
    elif isinstance(expr, (float, int)):
        return {"type": "Number", "value": expr}
    elif isinstance(expr, sympy.Float):
        return {"type": "Number", "value": float(expr)}
    elif isinstance(expr, sympy.Integer):
        return {"type": "Number", "value": int(expr)}
    else:
        raise NotImplementedError("Type not implemented: " + str(type(expr)))

if a[-1][1] is S.One:
        # there are no surds
        return p
    surds = [z for y, z in a]
    for i in range(len(surds)):
        if surds[i] != 1:
            break
    g, b1, b2 = _split_gcd(*surds[i:])
    a1 = []
    a2 = []
    for y, z in a:
        if z in b1:
            a1.append(y*z**S.Half)
        else:
            a2.append(y*z**S.Half)
    p1 = Add(*a1)
    p2 = Add(*a2)
    p = _mexpand(p1**2) - _mexpand(p2**2)
    return p

elif isinstance(expr, (Add, Mul, Pow, And, Or,
                                Eq, Ne, Lt, Gt, Le, Ge)):
            # we reconstruct the arithmetic expressions using the annotated
            # arguments
            args = expr.args

            # we treat the first element
            a = args[0]
            a_new = self._annotate(a, **settings)
            expr_new = a_new

            # then we treat the rest
            for a in args[1:]:
                a_new = self._annotate(a, **settings)
                if isinstance(expr, Add):
                    expr_new = Add(expr_new, a_new)
                elif isinstance(expr, Mul):
                    expr_new = Mul(expr_new, a_new)
                elif isinstance(expr, Pow):
                    expr_new = Pow(expr_new, a_new)
                elif isinstance(expr, And):
                    expr_new = And(expr_new, a_new)
                elif isinstance(expr, Or):
                    expr_new = Or(expr_new, a_new)
                elif isinstance(expr, Eq):
                    expr_new = Eq(expr_new, a_new)
                elif isinstance(expr, Ne):
                    expr_new = Ne(expr_new, a_new)
                elif isinstance(expr, Lt):
                    expr_new = Lt(expr_new, a_new)
                elif isinstance(expr, Le):
                    expr_new = Le(expr_new, a_new)

def convert_add(add):
    if add.ADD():
        lh = convert_add(add.additive(0))
        rh = convert_add(add.additive(1))
        return sympy.Add(lh, rh, evaluate=False)
    elif add.SUB():
        lh = convert_add(add.additive(0))
        rh = convert_add(add.additive(1))
        return sympy.Add(lh, -1 * rh, evaluate=False)
    else:
        return convert_mp(add.mp())

def _mass_visitor(node):
            if isinstance(node, Add):
                return Add(*map(lambda x: _mass_visitor(x), node.args))
            else:
                if self._find_obj_by_type(node, Dt):
                    return node
                return 0

# that case, we can avoid temporary symbols if we ensure the
                # correct substitution order.
                if subsdict[v] in subsdict:
                    # (x, y) -> (y, x),  we need a temporary variable
                    x = Dummy('x')
                    subslist.append((k, x))
                    final_subs.append((x, v))
                else:
                    # (x, y) -> (y, a),  x->y must be done last
                    # but before temporary variables are resolved
                    final_subs.insert(0, (k, v))
            else:
                subslist.append((k, v))
        subslist.extend(final_subs)
        new_terms.append(term.subs(subslist))
    return Add(*new_terms)

from random import randint
from sympy.logic.boolalg import Or

matchpy = import_module("matchpy")

if matchpy:
    from matchpy import Arity, Operation, CommutativeOperation, AssociativeOperation, OneIdentityOperation, CustomConstraint, Pattern, ReplacementRule, ManyToOneReplacer
    from matchpy.expressions.functions import op_iter, create_operation_expression, op_len
    from sympy.integrals.rubi.symbol import WC
    from matchpy import is_match, replace_all

    Operation.register(Integral)
    Operation.register(Pow)
    OneIdentityOperation.register(Pow)

    Operation.register(Add)
    OneIdentityOperation.register(Add)
    CommutativeOperation.register(Add)
    AssociativeOperation.register(Add)

    Operation.register(Mul)
    OneIdentityOperation.register(Mul)
    CommutativeOperation.register(Mul)
    AssociativeOperation.register(Mul)

    Operation.register(exp)
    Operation.register(log)
    Operation.register(gamma)
    Operation.register(uppergamma)
    Operation.register(fresnels)
    Operation.register(fresnelc)
    Operation.register(erf)