def test_equivalent_internal_lines_VT2conjT2_ambiguous_order_AT():
    # These diagrams invokes _determine_ambiguous() because the
    # dummies can not be ordered unambiguously by the key alone
    i, j, k, l, m, n = symbols("i j k l m n", below_fermi=True, cls=Dummy)
    a, b, c, d, e, f = symbols("a b c d e f", above_fermi=True, cls=Dummy)
    p1, p2, p3, p4 = symbols("p1 p2 p3 p4", above_fermi=True, cls=Dummy)
    h1, h2, h3, h4 = symbols("h1 h2 h3 h4", below_fermi=True, cls=Dummy)

    from sympy.utilities.iterables import variations

    # atv(abcd)att(abij)att(cdij)
    template = atv(p1, p2, p3, p4) * att(p1, p2, i, j) * att(p3, p4, i, j)
    permutator = variations([a, b, c, d], 4)
    base = template.subs(zip([p1, p2, p3, p4], permutator.next()))
    for permut in permutator:
        subslist = zip([p1, p2, p3, p4], permut)
        expr = template.subs(subslist)
        assert substitute_dummies(expr) == substitute_dummies(base)
    template = atv(p1, p2, p3, p4) * att(p1, p2, j, i) * att(p3, p4, i, j)
    permutator = variations([a, b, c, d], 4)
    base = template.subs(zip([p1, p2, p3, p4], permutator.next()))
    for permut in permutator:
        subslist = zip([p1, p2, p3, p4], permut)
        expr = template.subs(subslist)
        assert substitute_dummies(expr) == substitute_dummies(base)

def test_index_permutations_with_dummies():
    a, b, c, d = symbols('a b c d')
    p, q, r, s = symbols('p q r s', cls=Dummy)
    f, g = map(Function, 'fg')
    P = PermutationOperator

    # No dummy substitution necessary
    expr = f(a, b, p, q) - f(b, a, p, q)
    assert simplify_index_permutations(
        expr, [P(a, b)]) == P(a, b)*f(a, b, p, q)

    # Cases where dummy substitution is needed
    expected = P(a, b)*substitute_dummies(f(a, b, p, q))

    expr = f(a, b, p, q) - f(b, a, q, p)
    result = simplify_index_permutations(expr, [P(a, b)])
    assert expected == substitute_dummies(result)

    expr = f(a, b, q, p) - f(b, a, p, q)
    result = simplify_index_permutations(expr, [P(a, b)])
    assert expected == substitute_dummies(result)

    # A case where nothing can be done
    expr = f(a, b, q, p) - g(b, a, p, q)
    result = simplify_index_permutations(expr, [P(a, b)])
    assert expr == result

def test_equivalent_internal_lines_VT2conjT2_ambiguous_order():
    # These diagrams invokes _determine_ambiguous() because the
    # dummies can not be ordered unambiguously by the key alone
    i, j, k, l, m, n = symbols('i j k l m n', below_fermi=True, cls=Dummy)
    a, b, c, d, e, f = symbols('a b c d e f', above_fermi=True, cls=Dummy)
    p1, p2, p3, p4 = symbols('p1 p2 p3 p4', above_fermi=True, cls=Dummy)
    h1, h2, h3, h4 = symbols('h1 h2 h3 h4', below_fermi=True, cls=Dummy)

    from sympy.utilities.iterables import variations

    v = Function('v')
    t = Function('t')
    dums = _get_ordered_dummies

    # v(abcd)t(abij)t(cdij)
    template = v(p1, p2, p3, p4)*t(p1, p2, i, j)*t(p3, p4, i, j)
    permutator = variations([a, b, c, d], 4)
    base = template.subs(zip([p1, p2, p3, p4], next(permutator)))
    for permut in permutator:
        subslist = zip([p1, p2, p3, p4], permut)
        expr = template.subs(subslist)
        assert dums(base) != dums(expr)
        assert substitute_dummies(expr) == substitute_dummies(base)
    template = v(p1, p2, p3, p4)*t(p1, p2, j, i)*t(p3, p4, i, j)
    permutator = variations([a, b, c, d], 4)
    base = template.subs(zip([p1, p2, p3, p4], next(permutator)))
    for permut in permutator:
        subslist = zip([p1, p2, p3, p4], permut)
        expr = template.subs(subslist)
        assert dums(base) != dums(expr)
        assert substitute_dummies(expr) == substitute_dummies(base)

def test_dummy_order_inner_outer_lines_VT1T1T1T1_AT():
    ii, jj = symbols('i j', below_fermi=True)
    aa, bb = symbols('a b', above_fermi=True)
    k, l = symbols('k l', below_fermi=True, cls=Dummy)
    c, d = symbols('c d', above_fermi=True, cls=Dummy)

    # Coupled-Cluster T2 terms with V*T1*T1*T1*T1
    # non-equivalent substitutions (change of sign)
    exprs = [
        # permut t <=> swapping external lines
        atv(k, l, c, d)*att(c, ii)*att(d, jj)*att(aa, k)*att(bb, l),
        atv(k, l, c, d)*att(c, jj)*att(d, ii)*att(aa, k)*att(bb, l),
        atv(k, l, c, d)*att(c, ii)*att(d, jj)*att(bb, k)*att(aa, l),
    ]
    for permut in exprs[1:]:
        assert substitute_dummies(exprs[0]) == -substitute_dummies(permut)

    # equivalent substitutions
    exprs = [
        atv(k, l, c, d)*att(c, ii)*att(d, jj)*att(aa, k)*att(bb, l),
        # permut t <=> swapping external lines
        atv(k, l, c, d)*att(c, jj)*att(d, ii)*att(bb, k)*att(aa, l),
    ]
    for permut in exprs[1:]:
        assert substitute_dummies(exprs[0]) == substitute_dummies(permut)

def test_internal_external_pqrs_AT():
    ii, jj = symbols('i j')
    aa, bb = symbols('a b')
    k, l = symbols('k l', cls=Dummy)
    c, d = symbols('c d', cls=Dummy)

    exprs = [
        atv(k, l, c, d)*att(aa, c, ii, k)*att(bb, d, jj, l),
        atv(l, k, c, d)*att(aa, c, ii, l)*att(bb, d, jj, k),
        atv(k, l, d, c)*att(aa, d, ii, k)*att(bb, c, jj, l),
        atv(l, k, d, c)*att(aa, d, ii, l)*att(bb, c, jj, k),
    ]
    for permut in exprs[1:]:
        assert substitute_dummies(exprs[0]) == substitute_dummies(permut)

def test_equivalent_internal_lines_VT1T1_AT():
    i, j, k, l = symbols('i j k l', below_fermi=True, cls=Dummy)
    a, b, c, d = symbols('a b c d', above_fermi=True, cls=Dummy)

    exprs = [  # permute v.  Different dummy order. Not equivalent.
        atv(i, j, a, b)*att(a, i)*att(b, j),
        atv(j, i, a, b)*att(a, i)*att(b, j),
        atv(i, j, b, a)*att(a, i)*att(b, j),
    ]
    for permut in exprs[1:]:
        assert substitute_dummies(exprs[0]) != substitute_dummies(permut)

    exprs = [  # permute v.  Different dummy order. Equivalent
        atv(i, j, a, b)*att(a, i)*att(b, j),
        atv(j, i, b, a)*att(a, i)*att(b, j),
    ]
    for permut in exprs[1:]:
        assert substitute_dummies(exprs[0]) == substitute_dummies(permut)

    exprs = [  # permute t.  Same dummy order, not equivalent.
        atv(i, j, a, b)*att(a, i)*att(b, j),
        atv(i, j, a, b)*att(b, i)*att(a, j),
    ]
    for permut in exprs[1:]:
        assert substitute_dummies(exprs[0]) != substitute_dummies(permut)

    exprs = [  # permute v and t.  Different dummy order, equivalent
        atv(i, j, a, b)*att(a, i)*att(b, j),
        atv(j, i, a, b)*att(a, j)*att(b, i),
        atv(i, j, b, a)*att(b, i)*att(a, j),
        atv(j, i, b, a)*att(b, j)*att(a, i),
    ]
    for permut in exprs[1:]:
        assert substitute_dummies(exprs[0]) == substitute_dummies(permut)

def test_substitute_dummies_substitution_order():
    i, j, k, l = symbols("i j k l", below_fermi=True, cls=Dummy)
    f = Function("f")
    from sympy.utilities.iterables import variations

    for permut in variations([i, j, k, l], 4):
        assert substitute_dummies(f(*permut) - f(i, j, k, l)) == 0

def test_dummy_order_inner_outer_lines_VT1T1T1_AT():
    ii = symbols('i', below_fermi=True)
    aa = symbols('a', above_fermi=True)
    k, l = symbols('k l', below_fermi=True, cls=Dummy)
    c, d = symbols('c d', above_fermi=True, cls=Dummy)

    # Coupled-Cluster T1 terms with V*T1*T1*T1
    # t^{a}_{k} t^{c}_{i} t^{d}_{l} v^{lk}_{dc}
    exprs = [
        # permut v and t <=> swapping internal lines, equivalent
        # irrespective of symmetries in v
        atv(k, l, c, d)*att(c, ii)*att(d, l)*att(aa, k),
        atv(l, k, c, d)*att(c, ii)*att(d, k)*att(aa, l),
        atv(k, l, d, c)*att(d, ii)*att(c, l)*att(aa, k),
        atv(l, k, d, c)*att(d, ii)*att(c, k)*att(aa, l),
    ]
    for permut in exprs[1:]:
        assert substitute_dummies(exprs[0]) == substitute_dummies(permut)

def test_internal_external_pqrs():
    ii, jj = symbols("i j")
    aa, bb = symbols("a b")
    k, l = symbols("k l", cls=Dummy)
    c, d = symbols("c d", cls=Dummy)

    v = Function("v")
    t = Function("t")
    dums = _get_ordered_dummies

    exprs = [
        v(k, l, c, d) * t(aa, c, ii, k) * t(bb, d, jj, l),
        v(l, k, c, d) * t(aa, c, ii, l) * t(bb, d, jj, k),
        v(k, l, d, c) * t(aa, d, ii, k) * t(bb, c, jj, l),
        v(l, k, d, c) * t(aa, d, ii, l) * t(bb, c, jj, k),
    ]
    for permut in exprs[1:]:
        assert dums(exprs[0]) != dums(permut)
        assert substitute_dummies(exprs[0]) == substitute_dummies(permut)

def test_internal_external_pqrs():
    ii, jj = symbols('i j')
    aa, bb = symbols('a b')
    k, l = symbols('k l', cls=Dummy)
    c, d = symbols('c d', cls=Dummy)

    v = Function('v')
    t = Function('t')
    dums = _get_ordered_dummies

    exprs = [
        v(k, l, c, d)*t(aa, c, ii, k)*t(bb, d, jj, l),
        v(l, k, c, d)*t(aa, c, ii, l)*t(bb, d, jj, k),
        v(k, l, d, c)*t(aa, d, ii, k)*t(bb, c, jj, l),
        v(l, k, d, c)*t(aa, d, ii, l)*t(bb, c, jj, k),
    ]
    for permut in exprs[1:]:
        assert dums(exprs[0]) != dums(permut)
        assert substitute_dummies(exprs[0]) == substitute_dummies(permut)

def test_equivalent_internal_lines_VT2_AT():
    i, j, k, l = symbols("i j k l", below_fermi=True, cls=Dummy)
    a, b, c, d = symbols("a b c d", above_fermi=True, cls=Dummy)

    exprs = [
        # permute v. Same dummy order, not equivalent.
        atv(i, j, a, b) * att(a, b, i, j),
        atv(j, i, a, b) * att(a, b, i, j),
        atv(i, j, b, a) * att(a, b, i, j),
    ]
    for permut in exprs[1:]:
        assert substitute_dummies(exprs[0]) != substitute_dummies(permut)

    exprs = [
        # permute t.
        atv(i, j, a, b) * att(a, b, i, j),
        atv(i, j, a, b) * att(b, a, i, j),
        atv(i, j, a, b) * att(a, b, j, i),
    ]
    for permut in exprs[1:]:
        assert substitute_dummies(exprs[0]) != substitute_dummies(permut)

    exprs = [  # permute v and t.  Relabelling of dummies should be equivalent.
        atv(i, j, a, b) * att(a, b, i, j),
        atv(j, i, a, b) * att(a, b, j, i),
        atv(i, j, b, a) * att(b, a, i, j),
        atv(j, i, b, a) * att(b, a, j, i),
    ]
    for permut in exprs[1:]:
        assert substitute_dummies(exprs[0]) == substitute_dummies(permut)

def test_internal_external_VT2T2_AT():
    ii, jj = symbols("i j", below_fermi=True)
    aa, bb = symbols("a b", above_fermi=True)
    k, l = symbols("k l", below_fermi=True, cls=Dummy)
    c, d = symbols("c d", above_fermi=True, cls=Dummy)

    dums = _get_ordered_dummies

    exprs = [
        atv(k, l, c, d) * att(aa, c, ii, k) * att(bb, d, jj, l),
        atv(l, k, c, d) * att(aa, c, ii, l) * att(bb, d, jj, k),
        atv(k, l, d, c) * att(aa, d, ii, k) * att(bb, c, jj, l),
        atv(l, k, d, c) * att(aa, d, ii, l) * att(bb, c, jj, k),
    ]
    for permut in exprs[1:]:
        assert substitute_dummies(exprs[0]) == substitute_dummies(permut)
    exprs = [
        atv(k, l, c, d) * att(aa, c, ii, k) * att(d, bb, jj, l),
        atv(l, k, c, d) * att(aa, c, ii, l) * att(d, bb, jj, k),
        atv(k, l, d, c) * att(aa, d, ii, k) * att(c, bb, jj, l),
        atv(l, k, d, c) * att(aa, d, ii, l) * att(c, bb, jj, k),
    ]
    for permut in exprs[1:]:
        assert substitute_dummies(exprs[0]) == substitute_dummies(permut)
    exprs = [
        atv(k, l, c, d) * att(c, aa, ii, k) * att(bb, d, jj, l),
        atv(l, k, c, d) * att(c, aa, ii, l) * att(bb, d, jj, k),
        atv(k, l, d, c) * att(d, aa, ii, k) * att(bb, c, jj, l),
        atv(l, k, d, c) * att(d, aa, ii, l) * att(bb, c, jj, k),
    ]
    for permut in exprs[1:]:
        assert substitute_dummies(exprs[0]) == substitute_dummies(permut)

def test_internal_external_VT2T2_AT():
    ii, jj = symbols('i j', below_fermi=True)
    aa, bb = symbols('a b', above_fermi=True)
    k, l = symbols('k l', below_fermi=True, cls=Dummy)
    c, d = symbols('c d', above_fermi=True, cls=Dummy)

    exprs = [
        atv(k, l, c, d)*att(aa, c, ii, k)*att(bb, d, jj, l),
        atv(l, k, c, d)*att(aa, c, ii, l)*att(bb, d, jj, k),
        atv(k, l, d, c)*att(aa, d, ii, k)*att(bb, c, jj, l),
        atv(l, k, d, c)*att(aa, d, ii, l)*att(bb, c, jj, k),
    ]
    for permut in exprs[1:]:
        assert substitute_dummies(exprs[0]) == substitute_dummies(permut)
    exprs = [
        atv(k, l, c, d)*att(aa, c, ii, k)*att(d, bb, jj, l),
        atv(l, k, c, d)*att(aa, c, ii, l)*att(d, bb, jj, k),
        atv(k, l, d, c)*att(aa, d, ii, k)*att(c, bb, jj, l),
        atv(l, k, d, c)*att(aa, d, ii, l)*att(c, bb, jj, k),
    ]
    for permut in exprs[1:]:
        assert substitute_dummies(exprs[0]) == substitute_dummies(permut)
    exprs = [
        atv(k, l, c, d)*att(c, aa, ii, k)*att(bb, d, jj, l),
        atv(l, k, c, d)*att(c, aa, ii, l)*att(bb, d, jj, k),
        atv(k, l, d, c)*att(d, aa, ii, k)*att(bb, c, jj, l),
        atv(l, k, d, c)*att(d, aa, ii, l)*att(bb, c, jj, k),
    ]
    for permut in exprs[1:]:
        assert substitute_dummies(exprs[0]) == substitute_dummies(permut)

def test_equivalent_internal_lines_VT2conjT2():
    # this diagram requires special handling in TCE
    i, j, k, l, m, n = symbols('i j k l m n', below_fermi=True, cls=Dummy)
    a, b, c, d, e, f = symbols('a b c d e f', above_fermi=True, cls=Dummy)
    p1, p2, p3, p4 = symbols('p1 p2 p3 p4', above_fermi=True, cls=Dummy)
    h1, h2, h3, h4 = symbols('h1 h2 h3 h4', below_fermi=True, cls=Dummy)

    from sympy.utilities.iterables import variations

    v = Function('v')
    t = Function('t')
    dums = _get_ordered_dummies

    # v(abcd)t(abij)t(ijcd)
    template = v(p1, p2, p3, p4)*t(p1, p2, i, j)*t(i, j, p3, p4)
    permutator = variations([a, b, c, d], 4)
    base = template.subs(zip([p1, p2, p3, p4], next(permutator)))
    for permut in permutator:
        subslist = zip([p1, p2, p3, p4], permut)
        expr = template.subs(subslist)
        assert dums(base) != dums(expr)
        assert substitute_dummies(expr) == substitute_dummies(base)
    template = v(p1, p2, p3, p4)*t(p1, p2, j, i)*t(j, i, p3, p4)
    permutator = variations([a, b, c, d], 4)
    base = template.subs(zip([p1, p2, p3, p4], next(permutator)))
    for permut in permutator:
        subslist = zip([p1, p2, p3, p4], permut)
        expr = template.subs(subslist)
        assert dums(base) != dums(expr)
        assert substitute_dummies(expr) == substitute_dummies(base)

    # v(abcd)t(abij)t(jicd)
    template = v(p1, p2, p3, p4)*t(p1, p2, i, j)*t(j, i, p3, p4)
    permutator = variations([a, b, c, d], 4)
    base = template.subs(zip([p1, p2, p3, p4], next(permutator)))
    for permut in permutator:
        subslist = zip([p1, p2, p3, p4], permut)
        expr = template.subs(subslist)
        assert dums(base) != dums(expr)
        assert substitute_dummies(expr) == substitute_dummies(base)
    template = v(p1, p2, p3, p4)*t(p1, p2, j, i)*t(i, j, p3, p4)
    permutator = variations([a, b, c, d], 4)
    base = template.subs(zip([p1, p2, p3, p4], next(permutator)))
    for permut in permutator:
        subslist = zip([p1, p2, p3, p4], permut)
        expr = template.subs(subslist)
        assert dums(base) != dums(expr)
        assert substitute_dummies(expr) == substitute_dummies(base)

def test_dummy_order_inner_outer_lines_VT1T1T1():
    ii = symbols("i", below_fermi=True)
    aa = symbols("a", above_fermi=True)
    k, l = symbols("k l", below_fermi=True, cls=Dummy)
    c, d = symbols("c d", above_fermi=True, cls=Dummy)

    v = Function("v")
    t = Function("t")
    dums = _get_ordered_dummies

    # Coupled-Cluster T1 terms with V*T1*T1*T1
    # t^{a}_{k} t^{c}_{i} t^{d}_{l} v^{lk}_{dc}
    exprs = [
        # permut v and t <=> swapping internal lines, equivalent
        # irrespective of symmetries in v
        v(k, l, c, d) * t(c, ii) * t(d, l) * t(aa, k),
        v(l, k, c, d) * t(c, ii) * t(d, k) * t(aa, l),
        v(k, l, d, c) * t(d, ii) * t(c, l) * t(aa, k),
        v(l, k, d, c) * t(d, ii) * t(c, k) * t(aa, l),
    ]
    for permut in exprs[1:]:
        assert dums(exprs[0]) != dums(permut)
        assert substitute_dummies(exprs[0]) == substitute_dummies(permut)

def test_equivalent_internal_lines_VT1T1():
    i,j,k,l = symbols('i j k l',below_fermi=True, cls=Dummy)
    a,b,c,d = symbols('a b c d',above_fermi=True, cls=Dummy)

    f = Function('f')
    v = Function('v')
    t = Function('t')
    dums = _get_ordered_dummies

    exprs = [ # permute v.  Different dummy order. Not equivalent.
            v(i, j, a, b)*t(a, i)*t(b, j),
            v(j, i, a, b)*t(a, i)*t(b, j),
            v(i, j, b, a)*t(a, i)*t(b, j),
            ]
    for permut in exprs[1:]:
        assert dums(exprs[0]) != dums(permut)
        assert substitute_dummies(exprs[0]) != substitute_dummies(permut)

    exprs = [ # permute v.  Different dummy order. Equivalent
            v(i, j, a, b)*t(a, i)*t(b, j),
            v(j, i, b, a)*t(a, i)*t(b, j),
            ]
    for permut in exprs[1:]:
        assert dums(exprs[0]) != dums(permut)
        assert substitute_dummies(exprs[0]) == substitute_dummies(permut)

    exprs = [ # permute t.  Same dummy order, not equivalent.
            v(i, j, a, b)*t(a, i)*t(b, j),
            v(i, j, a, b)*t(b, i)*t(a, j),
            ]
    for permut in exprs[1:]:
        assert dums(exprs[0]) == dums(permut)
        assert substitute_dummies(exprs[0]) != substitute_dummies(permut)

    exprs = [ # permute v and t.  Different dummy order, equivalent
            v(i, j, a, b)*t(a, i)*t(b, j),
            v(j, i, a, b)*t(a, j)*t(b, i),
            v(i, j, b, a)*t(b, i)*t(a, j),
            v(j, i, b, a)*t(b, j)*t(a, i),
            ]
    for permut in exprs[1:]:
        assert dums(exprs[0]) != dums(permut)
        assert substitute_dummies(exprs[0]) == substitute_dummies(permut)

def test_dummy_order_inner_outer_lines_VT1T1T1T1():
    ii, jj = symbols('i j', below_fermi=True)
    aa, bb = symbols('a b', above_fermi=True)
    k, l = symbols('k l', below_fermi=True, cls=Dummy)
    c, d = symbols('c d', above_fermi=True, cls=Dummy)

    v = Function('v')
    t = Function('t')
    dums = _get_ordered_dummies

    # Coupled-Cluster T2 terms with V*T1*T1*T1*T1
    exprs = [
        # permut t <=> swapping external lines, not equivalent
        # except if v has certain symmetries.
        v(k, l, c, d)*t(c, ii)*t(d, jj)*t(aa, k)*t(bb, l),
        v(k, l, c, d)*t(c, jj)*t(d, ii)*t(aa, k)*t(bb, l),
        v(k, l, c, d)*t(c, ii)*t(d, jj)*t(bb, k)*t(aa, l),
        v(k, l, c, d)*t(c, jj)*t(d, ii)*t(bb, k)*t(aa, l),
    ]
    for permut in exprs[1:]:
        assert dums(exprs[0]) != dums(permut)
        assert substitute_dummies(exprs[0]) != substitute_dummies(permut)
    exprs = [
        # permut v <=> swapping external lines, not equivalent
        # except if v has certain symmetries.
        #
        # Note that in contrast to above, these permutations have identical
        # dummy order.  That is because the proximity to external indices
        # has higher influence on the canonical dummy ordering than the
        # position of a dummy on the factors.  In fact, the terms here are
        # similar in structure as the result of the dummy substitions above.
        v(k, l, c, d)*t(c, ii)*t(d, jj)*t(aa, k)*t(bb, l),
        v(l, k, c, d)*t(c, ii)*t(d, jj)*t(aa, k)*t(bb, l),
        v(k, l, d, c)*t(c, ii)*t(d, jj)*t(aa, k)*t(bb, l),
        v(l, k, d, c)*t(c, ii)*t(d, jj)*t(aa, k)*t(bb, l),
    ]
    for permut in exprs[1:]:
        assert dums(exprs[0]) == dums(permut)
        assert substitute_dummies(exprs[0]) != substitute_dummies(permut)
    exprs = [
        # permut t and v <=> swapping internal lines, equivalent.
        # Canonical dummy order is different, and a consistent
        # substitution reveals the equivalence.
        v(k, l, c, d)*t(c, ii)*t(d, jj)*t(aa, k)*t(bb, l),
        v(k, l, d, c)*t(c, jj)*t(d, ii)*t(aa, k)*t(bb, l),
        v(l, k, c, d)*t(c, ii)*t(d, jj)*t(bb, k)*t(aa, l),
        v(l, k, d, c)*t(c, jj)*t(d, ii)*t(bb, k)*t(aa, l),
    ]
    for permut in exprs[1:]:
        assert dums(exprs[0]) != dums(permut)
        assert substitute_dummies(exprs[0]) == substitute_dummies(permut)

def test_equivalent_internal_lines_VT2conjT2_AT():
    # this diagram requires special handling in TCE
    i, j, k, l, m, n = symbols("i j k l m n", below_fermi=True, cls=Dummy)
    a, b, c, d, e, f = symbols("a b c d e f", above_fermi=True, cls=Dummy)
    p1, p2, p3, p4 = symbols("p1 p2 p3 p4", above_fermi=True, cls=Dummy)
    h1, h2, h3, h4 = symbols("h1 h2 h3 h4", below_fermi=True, cls=Dummy)

    from sympy.utilities.iterables import variations

    # atv(abcd)att(abij)att(ijcd)
    template = atv(p1, p2, p3, p4) * att(p1, p2, i, j) * att(i, j, p3, p4)
    permutator = variations([a, b, c, d], 4)
    base = template.subs(zip([p1, p2, p3, p4], permutator.next()))
    for permut in permutator:
        subslist = zip([p1, p2, p3, p4], permut)
        expr = template.subs(subslist)
        assert substitute_dummies(expr) == substitute_dummies(base)
    template = atv(p1, p2, p3, p4) * att(p1, p2, j, i) * att(j, i, p3, p4)
    permutator = variations([a, b, c, d], 4)
    base = template.subs(zip([p1, p2, p3, p4], permutator.next()))
    for permut in permutator:
        subslist = zip([p1, p2, p3, p4], permut)
        expr = template.subs(subslist)
        assert substitute_dummies(expr) == substitute_dummies(base)

    # atv(abcd)att(abij)att(jicd)
    template = atv(p1, p2, p3, p4) * att(p1, p2, i, j) * att(j, i, p3, p4)
    permutator = variations([a, b, c, d], 4)
    base = template.subs(zip([p1, p2, p3, p4], permutator.next()))
    for permut in permutator:
        subslist = zip([p1, p2, p3, p4], permut)
        expr = template.subs(subslist)
        assert substitute_dummies(expr) == substitute_dummies(base)
    template = atv(p1, p2, p3, p4) * att(p1, p2, j, i) * att(i, j, p3, p4)
    permutator = variations([a, b, c, d], 4)
    base = template.subs(zip([p1, p2, p3, p4], permutator.next()))
    for permut in permutator:
        subslist = zip([p1, p2, p3, p4], permut)
        expr = template.subs(subslist)
        assert substitute_dummies(expr) == substitute_dummies(base)

def test_equivalent_internal_lines_VT2():
    i,j,k,l = symbols('i j k l',below_fermi=True, cls=Dummy)
    a,b,c,d = symbols('a b c d',above_fermi=True, cls=Dummy)

    f = Function('f')
    v = Function('v')
    t = Function('t')
    dums = _get_ordered_dummies
    exprs = [
            # permute v. Same dummy order, not equivalent.
            #
            # This test show that the dummy order may not be sensitive to all
            # index permutations.  The following expressions have identical
            # structure as the resulting terms from of the dummy subsitutions
            # in the test above.  Here, all expressions have the same dummy
            # order, so they cannot be simplified by means of dummy
            # substitution.  In order to simplify further, it is necessary to
            # exploit symmetries in the objects, for instance if t or v is
            # antisymmetric.
            v(i, j, a, b)*t(a, b, i, j),
            v(j, i, a, b)*t(a, b, i, j),
            v(i, j, b, a)*t(a, b, i, j),
            v(j, i, b, a)*t(a, b, i, j),
            ]
    for permut in exprs[1:]:
        assert dums(exprs[0]) == dums(permut)
        assert substitute_dummies(exprs[0]) != substitute_dummies(permut)

    exprs = [
            # permute t.
            v(i, j, a, b)*t(a, b, i, j),
            v(i, j, a, b)*t(b, a, i, j),
            v(i, j, a, b)*t(a, b, j, i),
            v(i, j, a, b)*t(b, a, j, i),
            ]
    for permut in exprs[1:]:
        assert dums(exprs[0]) != dums(permut)
        assert substitute_dummies(exprs[0]) != substitute_dummies(permut)

    exprs = [ # permute v and t.  Relabelling of dummies should be equivalent.
            v(i, j, a, b)*t(a, b, i, j),
            v(j, i, a, b)*t(a, b, j, i),
            v(i, j, b, a)*t(b, a, i, j),
            v(j, i, b, a)*t(b, a, j, i),
            ]
    for permut in exprs[1:]:
        assert dums(exprs[0]) != dums(permut)
        assert substitute_dummies(exprs[0]) == substitute_dummies(permut)

def test_dummy_order_ambiguous():
    aa, bb = symbols('ab', above_fermi=True)
    i, j, k, l, m = symbols('i j k l m', below_fermi=True, cls=Dummy)
    a, b, c, d, e = symbols('a b c d e', above_fermi=True, cls=Dummy)
    p, q = symbols('p q', cls=Dummy)
    p1,p2,p3,p4 = symbols('p1 p2 p3 p4',above_fermi=True, cls=Dummy)
    p5,p6,p7,p8 = symbols('p5 p6 p7 p8',above_fermi=True, cls=Dummy)
    h1,h2,h3,h4 = symbols('h1 h2 h3 h4',below_fermi=True, cls=Dummy)
    h5,h6,h7,h8 = symbols('h5 h6 h7 h8',below_fermi=True, cls=Dummy)

    A = Function('A')
    B = Function('B')
    dums = _get_ordered_dummies

    from sympy.utilities.iterables import variations

    # A*A*A*A*B  --  ordering of p5 and p4 is used to figure out the rest
    template = A(p1, p2)*A(p4, p1)*A(p2, p3)*A(p3, p5)*B(p5, p4)
    permutator = variations([a,b,c,d,e], 5)
    base = template.subs(zip([p1, p2, p3, p4, p5], permutator.next()))
    for permut in permutator:
        subslist = zip([p1, p2, p3, p4, p5], permut)
        expr = template.subs(subslist)
        assert substitute_dummies(expr) == substitute_dummies(base)

    # A*A*A*A*A  --  an arbitrary index is assigned and the rest are figured out
    template = A(p1, p2)*A(p4, p1)*A(p2, p3)*A(p3, p5)*A(p5, p4)
    permutator = variations([a,b,c,d,e], 5)
    base = template.subs(zip([p1, p2, p3, p4, p5], permutator.next()))
    for permut in permutator:
        subslist = zip([p1, p2, p3, p4, p5], permut)
        expr = template.subs(subslist)
        assert substitute_dummies(expr) == substitute_dummies(base)

    # A*A*A  --  ordering of p5 and p4 is used to figure out the rest
    template = A(p1, p2, p4, p1)*A(p2, p3, p3, p5)*A(p5, p4)
    permutator = variations([a,b,c,d,e], 5)
    base = template.subs(zip([p1, p2, p3, p4, p5], permutator.next()))
    for permut in permutator:
        subslist = zip([p1, p2, p3, p4, p5], permut)
        expr = template.subs(subslist)
        assert substitute_dummies(expr) == substitute_dummies(base)