Backport change in reposition_dummy_indices to 1.3.

diff --git a/check/exam_paranoia.cpp b/check/exam_paranoia.cpp
index 7c82e48662d9181256a5ebc044edd2b501f3d3e4..0ab64910c437e7d8c51d773cdacc5637680307b5 100644 (file)
--- a/check/exam_paranoia.cpp
+++ b/check/exam_paranoia.cpp
@@ -437,6 +437,33 @@ static unsigned exam_paranoia16()
        return result;
 }
 
+// Bug in reposition_dummy_indices() could result in correct expression
+// turned into one with inconsistent indices. Fixed on Aug 29, 2006
+static unsigned exam_paranoia17()
+{
+       varidx mu1(symbol("mu1"), 4);
+       varidx mu2(symbol("mu2"), 4);
+       varidx mu3(symbol("mu3"), 4);
+       varidx mu4(symbol("mu4"), 4);
+       varidx mu5(symbol("mu5"), 4);
+       varidx mu6(symbol("mu6"), 4);
+
+       exvector ev2;
+       ev2.push_back(mu3.toggle_variance());
+       ev2.push_back(mu6);
+       ev2.push_back(mu5.toggle_variance());
+       ev2.push_back(mu6.toggle_variance());
+       ev2.push_back(mu5);
+       ev2.push_back(mu3); 
+       // notice: all indices are contracted ...
+
+       ex test_cycl = indexed(symbol("A"), sy_cycl(), ev2);
+       test_cycl = test_cycl.simplify_indexed();
+       // ... so there should be zero free indices in the end.
+       return test_cycl.get_free_indices().size();
+}
+
+
 unsigned exam_paranoia()
 {
        unsigned result = 0;
@@ -460,6 +487,7 @@ unsigned exam_paranoia()
        result += exam_paranoia14();  cout << '.' << flush;
        result += exam_paranoia15();  cout << '.' << flush;
        result += exam_paranoia16();  cout << '.' << flush;
+       result += exam_paranoia17();  cout << '.' << flush;
        
        if (!result) {
                cout << " passed " << endl;

diff --git a/ginac/indexed.cpp b/ginac/indexed.cpp
index 05c3b90a58f0fc793344d4a02c7eb69718957edb..2de726bfdb2e46bf612c3ce29696493f19209b24 100644 (file)
--- a/ginac/indexed.cpp
+++ b/ginac/indexed.cpp
@@ -632,10 +632,60 @@ bool reposition_dummy_indices(ex & e, exvector & variant_dummy_indices, exvector
 {
        bool something_changed = false;
 
+       // Find dummy symbols that occur twice in the same indexed object.
+       exvector local_var_dummies;
+       local_var_dummies.reserve(e.nops()/2);
+       for (size_t i=1; i<e.nops(); ++i) {
+               if (!is_a<varidx>(e.op(i)))
+                       continue;
+               for (size_t j=i+1; j<e.nops(); ++j) {
+                       if (is_dummy_pair(e.op(i), e.op(j))) {
+                               local_var_dummies.push_back(e.op(i));
+                               for (exvector::iterator k = variant_dummy_indices.begin();
+                                               k!=variant_dummy_indices.end(); ++k) {
+                                       if (e.op(i).op(0) == k->op(0)) {
+                                               variant_dummy_indices.erase(k);
+                                               break;
+                                       }
+                               }
+                               break;
+                       }
+               }
+       }
+
+       // In the case where a dummy symbol occurs twice in the same indexed object
+       // we try all posibilities of raising/lowering and keep the least one in
+       // the sense of ex_is_less.
+       ex optimal_e = e;
+       size_t numpossibs = 1 << local_var_dummies.size();
+       for (size_t i=0; i<numpossibs; ++i) {
+               ex try_e = e;
+               for (size_t j=0; j<local_var_dummies.size(); ++j) {
+                       exmap m;
+                       if (1<<j & i) {
+                               ex curr_idx = local_var_dummies[j];
+                               ex curr_toggle = ex_to<varidx>(curr_idx).toggle_variance();
+                               m[curr_idx] = curr_toggle;
+                               m[curr_toggle] = curr_idx;
+                       }
+                       try_e = e.subs(m, subs_options::no_pattern);
+               }
+               if(ex_is_less()(try_e, optimal_e))
+               {       optimal_e = try_e;
+                       something_changed = true;
+               }
+       }
+       e = optimal_e;
+
+       if (!is_a<indexed>(e))
+               return true;
+
+       exvector seq = ex_to<indexed>(e).seq;
+
        // If a dummy index is encountered for the first time in the
        // product, pull it up, otherwise, pull it down
-       exvector::const_iterator it2, it2start, it2end;
-       for (it2start = ex_to<indexed>(e).seq.begin(), it2end = ex_to<indexed>(e).seq.end(), it2 = it2start + 1; it2 != it2end; ++it2) {
+       for (exvector::iterator it2 = seq.begin()+1, it2end = seq.end();
+                       it2 != it2end; ++it2) {
                if (!is_exactly_a<varidx>(*it2))
                        continue;
 
@@ -643,14 +693,20 @@ bool reposition_dummy_indices(ex & e, exvector & variant_dummy_indices, exvector
                for (vit = variant_dummy_indices.begin(), vitend = variant_dummy_indices.end(); vit != vitend; ++vit) {
                        if (it2->op(0).is_equal(vit->op(0))) {
                                if (ex_to<varidx>(*it2).is_covariant()) {
-                                       e = e.subs(lst(
-                                               *it2 == ex_to<varidx>(*it2).toggle_variance(),
-                                               ex_to<varidx>(*it2).toggle_variance() == *it2
-                                       ), subs_options::no_pattern);
+                                       /*
+                                        * N.B. we don't want to use
+                                        *
+                                        *  e = e.subs(lst(
+                                        *  *it2 == ex_to<varidx>(*it2).toggle_variance(),
+                                        *  ex_to<varidx>(*it2).toggle_variance() == *it2
+                                        *  ), subs_options::no_pattern);
+                                        *
+                                        * since this can trigger non-trivial repositioning of indices,
+                                        * e.g. due to non-trivial symmetry properties of e, thus
+                                        * invalidating iterators
+                                        */
+                                       *it2 = ex_to<varidx>(*it2).toggle_variance();
                                        something_changed = true;
-                                       it2 = ex_to<indexed>(e).seq.begin() + (it2 - it2start);
-                                       it2start = ex_to<indexed>(e).seq.begin();
-                                       it2end = ex_to<indexed>(e).seq.end();
                                }
                                moved_indices.push_back(*vit);
                                variant_dummy_indices.erase(vit);
@@ -661,11 +717,8 @@ bool reposition_dummy_indices(ex & e, exvector & variant_dummy_indices, exvector
                for (vit = moved_indices.begin(), vitend = moved_indices.end(); vit != vitend; ++vit) {
                        if (it2->op(0).is_equal(vit->op(0))) {
                                if (ex_to<varidx>(*it2).is_contravariant()) {
-                                       e = e.subs(*it2 == ex_to<varidx>(*it2).toggle_variance(), subs_options::no_pattern);
+                                       *it2 = ex_to<varidx>(*it2).toggle_variance();
                                        something_changed = true;
-                                       it2 = ex_to<indexed>(e).seq.begin() + (it2 - it2start);
-                                       it2start = ex_to<indexed>(e).seq.begin();
-                                       it2end = ex_to<indexed>(e).seq.end();
                                }
                                goto next_index;
                        }
@@ -674,6 +727,9 @@ bool reposition_dummy_indices(ex & e, exvector & variant_dummy_indices, exvector
 next_index: ;
        }
 
+       if (something_changed)
+               e = ex_to<indexed>(e).thiscontainer(seq);
+
        return something_changed;
 }