Execute sums over KroneckerDelta prefactors

def indexed_sum_over_kronecker(cls, ops, kwargs):
    """Execute sums over KroneckerDelta prefactors"""
    from qnet.algebra.core.abstract_quantum_algebra import QuantumExpression
    term, *ranges = ops
    assert isinstance(term, QuantumExpression)
    deltas = set(Pattern(head=sympy.KroneckerDelta).findall(term))
    if len(deltas) == 0:
        return ops, kwargs  # nothing to do
    else:  # the term contains at least one KroneckerDelta
        delta_indices = set.union(*[set(
            [idx for idx in delta.free_symbols if isinstance(idx, IdxSym)])
            for delta in deltas])
        ranges = sorted(  # sort in the order we'd prefer to eliminate
            ranges,
            key=partial(_ranges_key, delta_indices=delta_indices),
            reverse=True)
        buffer = [(term, ranges)]
        i = 0  # position in buffer that we're currently handling
        i_range = 0  # position of index-range for current buffer item
        while i < len(buffer):
            t, rs = buffer[i]
            if rs[i_range].index_symbol in delta_indices:
                new_items, flag = _deltasummation(t, rs, i_range)
                new_rs = new_items[0][1]  # same for all new_items
                buffer = buffer[:i] + new_items + buffer[i + 1:]
                assert flag in [1, 2, 3]
                if flag == 2:
                    i_range += 1
                    # * for flag == 1, leaving i_range unchanged will
                    # effectively to to the next range (as the current range
                    # was removed)
                    # * for flag == 3, buffer[i] has changed, and we'll want to
                    # call it again with the same i_range
            else:
                # if the index symbol doesn't occur in any KroneckerDelta,
                # there is no chance _deltasummation will do anything; so we
                # just skip to the next index
                i_range += 1
                new_rs = rs
            if i_range >= len(new_rs):
                # if we've exhausted the index-ranges for the current buffer
                # item, go to the next buffer item
                i += 1
                i_range = 0
        if len(buffer) == 1 and buffer[0] == (term, ranges):
            return ops, kwargs  # couldn't resolve deltas
        else:
            (t, rs) = buffer[0]
            res = t
            if len(rs) > 0:
                res = cls.create(t, *rs, **kwargs)
            for (t, rs) in buffer[1:]:
                if len(rs) > 0:
                    t = cls.create(t, *rs, **kwargs)
                res += t
            return res