Index: ginac/mul.cpp =================================================================== RCS file: /home/cvs/GiNaC/ginac/mul.cpp,v retrieving revision 1.90 diff -u -r1.90 mul.cpp --- ginac/mul.cpp 24 Jul 2005 21:02:39 -0000 1.90 +++ ginac/mul.cpp 16 Nov 2005 21:50:30 -0000 @@ -895,16 +895,33 @@ // Compute the new overall coefficient and put it together: ex tmp_accu = (new add(distrseq, add1.overall_coeff*add2.overall_coeff))->setflag(status_flags::dynallocated); + exvector add1_dummy_indices, add2_dummy_indices, add_indices; + + for (epvector::const_iterator i=add1begin; i!=add1end; ++i) { + add_indices = get_all_dummy_indices(i->rest); + add1_dummy_indices.insert(add1_dummy_indices.end(), add_indices.begin(), add_indices.end()); + } + for (epvector::const_iterator i=add2begin; i!=add2end; ++i) { + add_indices = get_all_dummy_indices(i->rest); + add2_dummy_indices.insert(add2_dummy_indices.end(), add_indices.begin(), add_indices.end()); + } + + sort(add1_dummy_indices.begin(), add1_dummy_indices.end(), ex_is_less()); + sort(add2_dummy_indices.begin(), add2_dummy_indices.end(), ex_is_less()); + lst dummy_subs = rename_dummy_indices_uniquely(add1_dummy_indices, add2_dummy_indices); + // Multiply explicitly all non-numeric terms of add1 and add2: - for (epvector::const_iterator i1=add1begin; i1!=add1end; ++i1) { + for (epvector::const_iterator i2=add2begin; i2!=add2end; ++i2) { // We really have to combine terms here in order to compactify // the result. Otherwise it would become waayy tooo bigg. numeric oc; distrseq.clear(); - for (epvector::const_iterator i2=add2begin; i2!=add2end; ++i2) { + ex i2_new = (dummy_subs.op(0).nops()>0? + i2->rest.subs((lst)dummy_subs.op(0), (lst)dummy_subs.op(1), subs_options::no_pattern) : i2->rest); + for (epvector::const_iterator i1=add1begin; i1!=add1end; ++i1) { // Don't push_back expairs which might have a rest that evaluates to a numeric, // since that would violate an invariant of expairseq: - const ex rest = (new mul(i1->rest, rename_dummy_indices_uniquely(i1->rest, i2->rest)))->setflag(status_flags::dynallocated); + const ex rest = (new mul(i1->rest, i2_new))->setflag(status_flags::dynallocated); if (is_exactly_a(rest)) { oc += ex_to(rest).mul(ex_to(i1->coeff).mul(ex_to(i2->coeff))); } else { @@ -932,10 +949,12 @@ size_t n = last_expanded.nops(); exvector distrseq; distrseq.reserve(n); + exvector va = get_all_dummy_indices(mul(non_adds)); + sort(va.begin(), va.end(), ex_is_less()); for (size_t i=0; isetflag(status_flags::dynallocated); if (can_be_further_expanded(term)) { distrseq.push_back(term.expand()); Index: ginac/ncmul.cpp =================================================================== RCS file: /home/cvs/GiNaC/ginac/ncmul.cpp,v retrieving revision 1.57 diff -u -r1.57 ncmul.cpp --- ginac/ncmul.cpp 19 May 2005 14:10:40 -0000 1.57 +++ ginac/ncmul.cpp 16 Nov 2005 21:50:30 -0000 @@ -170,21 +170,22 @@ /* Rename indices in the static members of the product */ exvector expanded_seq_mod; - size_t j = 0; - size_t i; + size_t j = 0, i; + exvector va; + for (i=0; i Index: ginac/power.cpp =================================================================== RCS file: /home/cvs/GiNaC/ginac/power.cpp,v retrieving revision 1.102 diff -u -r1.102 power.cpp --- ginac/power.cpp 9 Nov 2005 13:47:00 -0000 1.102 +++ ginac/power.cpp 16 Nov 2005 21:50:30 -0000 @@ -864,8 +864,11 @@ // Leave it to multiplication since dummy indices have to be renamed if (get_all_dummy_indices(m).size() > 0 && n.is_positive()) { ex result = m; + exvector va = get_all_dummy_indices(m); + sort(va.begin(), va.end(), ex_is_less()); + for (int i=1; i < n.to_int(); i++) - result *= rename_dummy_indices_uniquely(m,m); + result *= rename_dummy_indices_uniquely(va, m); return result; } Index: ginac/indexed.h =================================================================== RCS file: /home/cvs/GiNaC/ginac/indexed.h,v retrieving revision 1.53 diff -u -r1.53 indexed.h --- ginac/indexed.h 11 Nov 2005 12:05:24 -0000 1.53 +++ ginac/indexed.h 16 Nov 2005 21:50:30 -0000 @@ -254,12 +254,20 @@ /** Returns all dummy indices from the expression */ exvector get_all_dummy_indices(const ex & e); +/** Returns b with all dummy indices, which are listed in va, renamed + * if modify_va is set to TRUE all dummy indices of b will be appended to va */ +ex rename_dummy_indices_uniquely(exvector & va, const ex & b, bool modify_va = false); + /** Returns b with all dummy indices, which are common with a, renamed */ ex rename_dummy_indices_uniquely(const ex & a, const ex & b); /** Same as above, where va and vb contain the indices of a and b and are sorted */ ex rename_dummy_indices_uniquely(const exvector & va, const exvector & vb, const ex & b); +/** Similar to above, where va and vb are the same and the return value is a list of two lists + * for substitution in b */ +lst rename_dummy_indices_uniquely(const exvector & va, const exvector & vb); + /** This function returns the given expression with expanded sums * for all dummy index summations, where the dimensionality of * the dummy index is numeric. Index: ginac/indexed.cpp =================================================================== RCS file: /home/cvs/GiNaC/ginac/indexed.cpp,v retrieving revision 1.97 diff -u -r1.97 indexed.cpp --- ginac/indexed.cpp 11 Nov 2005 12:05:24 -0000 1.97 +++ ginac/indexed.cpp 16 Nov 2005 21:50:30 -0000 @@ -1374,12 +1374,12 @@ return v; } -ex rename_dummy_indices_uniquely(const exvector & va, const exvector & vb, const ex & b) +lst rename_dummy_indices_uniquely(const exvector & va, const exvector & vb) { exvector common_indices; set_intersection(va.begin(), va.end(), vb.begin(), vb.end(), std::back_insert_iterator(common_indices), ex_is_less()); if (common_indices.empty()) { - return b; + return lst(lst(), lst()); } else { exvector new_indices, old_indices; old_indices.reserve(2*common_indices.size()); @@ -1408,17 +1408,57 @@ } ++ip; } - return b.subs(lst(old_indices.begin(), old_indices.end()), lst(new_indices.begin(), new_indices.end()), subs_options::no_pattern); + return lst(lst(old_indices.begin(), old_indices.end()), lst(new_indices.begin(), new_indices.end())); } } +ex rename_dummy_indices_uniquely(const exvector & va, const exvector & vb, const ex & b) +{ + lst indices_subs = rename_dummy_indices_uniquely(va, vb); + return (indices_subs.op(0).nops()>0 ? b.subs((lst)indices_subs.op(0), (lst)indices_subs.op(1), subs_options::no_pattern) : b); +} + ex rename_dummy_indices_uniquely(const ex & a, const ex & b) { exvector va = get_all_dummy_indices(a); - exvector vb = get_all_dummy_indices(b); - sort(va.begin(), va.end(), ex_is_less()); - sort(vb.begin(), vb.end(), ex_is_less()); - return rename_dummy_indices_uniquely(va, vb, b); + if (va.size() > 0) { + exvector vb = get_all_dummy_indices(b); + if (vb.size() > 0) { + sort(va.begin(), va.end(), ex_is_less()); + sort(vb.begin(), vb.end(), ex_is_less()); + lst indices_subs = rename_dummy_indices_uniquely(va, vb); + if (indices_subs.op(0).nops() > 0) + return b.subs((lst)indices_subs.op(0), (lst)indices_subs.op(1), subs_options::no_pattern); + } + } + return b; +} + +ex rename_dummy_indices_uniquely(exvector & va, const ex & b, bool modify_va) +{ + if (va.size() > 0) { + exvector vb = get_all_dummy_indices(b); + if (vb.size() > 0) { + sort(vb.begin(), vb.end(), ex_is_less()); + lst indices_subs = rename_dummy_indices_uniquely(va, vb); + if (indices_subs.op(0).nops() > 0) { + if (modify_va) { + for (lst::const_iterator i = ex_to(indices_subs.op(1)).begin(); i != ex_to(indices_subs.op(1)).end(); ++i) + va.push_back(*i); + exvector uncommon_indices; + set_difference(vb.begin(), vb.end(), indices_subs.op(0).begin(), indices_subs.op(0).end(), std::back_insert_iterator(uncommon_indices), ex_is_less()); + exvector::const_iterator ip = uncommon_indices.begin(), ipend = uncommon_indices.end(); + while (ip != ipend) { + va.push_back(*ip); + ++ip; + } + sort(va.begin(), va.end(), ex_is_less()); + } + return b.subs((lst)indices_subs.op(0), (lst)indices_subs.op(1), subs_options::no_pattern); + } + } + } + return b; } ex expand_dummy_sum(const ex & e, bool subs_idx) Index: check/exam_clifford.cpp =================================================================== RCS file: /home/cvs/GiNaC/check/exam_clifford.cpp,v retrieving revision 1.27 diff -u -r1.27 exam_clifford.cpp --- check/exam_clifford.cpp 12 Jul 2005 17:56:29 -0000 1.27 +++ check/exam_clifford.cpp 16 Nov 2005 21:50:30 -0000 @@ -48,9 +48,9 @@ static unsigned check_equal_lst(const ex & e1, const ex & e2) { - for(int i = 0; i++; i < e1.nops()) { + for (unsigned int i = 0; i < e1.nops(); i++) { ex e = e1.op(i) - e2.op(i); - if (!e.is_zero()) { + if (!e.normal().is_zero()) { clog << "(" << e1 << ") - (" << e2 << ") erroneously returned " << e << " instead of 0 (in the entry " << i << ")" << endl; return 1; @@ -314,7 +314,7 @@ matrix A_symm(4,4), A2(4, 4); A_symm = A.add(A.transpose()).mul(half); A2 = A_symm.mul(A_symm); - + ex e, e1; bool anticommuting = ex_to(clifford_unit(nu, A)).is_anticommuting(); int result = 0; @@ -412,7 +412,7 @@ ex c = clifford_unit(nu, A, 1); e = lst_to_clifford(lst(t, x, y, z), mu, A, 1) * lst_to_clifford(lst(1, 2, 3, 4), c); e1 = clifford_inverse(e); - result += check_equal_lst((e*e1).simplify_indexed(), dirac_ONE(1)); + result += check_equal((e*e1).simplify_indexed(), dirac_ONE(1)); // Moebius map (both forms) checks for symmetric metrics only matrix M1(2, 2), M2(2, 2);