Expand the first :class:`sympy.Add` containing a simple
    :class:`sympy.KroneckerDelta`.

    Auxiliary routine for :func:`_deltasummation`. Adapted from SymPy. The
    input `expr` may be a :class:`.QuantumExpression` or a
    `:class:`sympy.Basic` instance.

    Returns a list of summands. The elements of the list may be
    :class:`.QuantumExpression` or a `:class:`sympy.Basic` instances. There is
    no guarantee of type stability: an input :class:`.QuantumExpression` may
    result in a :class:`sympy.Basic` instance in the `summands`.

def _expand_delta(expr, idx):
    """Expand the first :class:`sympy.Add` containing a simple
    :class:`sympy.KroneckerDelta`.

    Auxiliary routine for :func:`_deltasummation`. Adapted from SymPy. The
    input `expr` may be a :class:`.QuantumExpression` or a
    `:class:`sympy.Basic` instance.

    Returns a list of summands. The elements of the list may be
    :class:`.QuantumExpression` or a `:class:`sympy.Basic` instances. There is
    no guarantee of type stability: an input :class:`.QuantumExpression` may
    result in a :class:`sympy.Basic` instance in the `summands`.
    """
    found_first_delta = False
    summands = None
    for factor in _factors_for_expand_delta(expr):
        need_to_expand = False
        if not found_first_delta and isinstance(factor, sympy.Basic):
            if factor.is_Add and _has_simple_delta(factor, idx):
                need_to_expand = True
        if need_to_expand:
            found_first_delta = True
            if summands is None:
                summands = list(factor.args)
            else:
                summands = [summands[0]*t for t in factor.args]
        else:
            if summands is None:
                summands = [factor, ]
            else:
                summands = [t*factor for t in summands]
    return summands