Here is a solution that solves some of these issues:
{a,b,c,d,a#}::Indices(vector).
{i,j,k,l,i#}::Indices(isospin).
\nabla{#}::Derivative.
\partial{#}::PartialDerivative.
def expand_nabla(ex):
    for nabla in ex[r'\nabla']:
        nabla.name=r'\partial'
        dindex = nabla.indices().__next__() 
        for arg in nabla.args():
            ret:=0;
            for index in arg.free_indices():
                indexprop = Indices.get(index, ignore_parent_rel = True)
                if indexprop.set_name == 'vector':
                    t2:= @(arg);
                    dummy = indexprop.get_dummy(ex)
                    if index.parent_rel==sub:
                        t1:= -\Gamma^{@(dummy)}_{@(dindex) @(index)};
                        t2[index]:= _{@(dummy)};
                    else:
                        t1:=  \Gamma^{@(index)}_{@(dindex) @(dummy)};
                        t2[index]:= ^{@(dummy)};
                    ret += Ex(str(nabla.multiplier)) * t1 * t2
            nabla += ret
    return ex
This version grabs a dummy index by looking at the entire expression tree, while in practice, we might just want to look at the term in the top-level sum. It knows to ignore isospin indices and only grab vector indices. We also need \nabla and \partial to have IndexInherit (inherited from Derivative).